PyIntact manual *************** Overview ======== * Introduction to the API This document provides a comprehensive guide to using our Python API for Intact.Simulation. The API allows users to define geometries, material properties, components, assemblies, and boundary conditions. Users can specify units, create simulation scenarios, execute the simulation, and query the results. The simulations can be run with different solver types and the results can be written to a file for visualization. The document also provides Python code examples and snippets for each step of the simulation process. * Key features and capabilities * Physics Support * Stress * Modal * Thermal * Geometry Support * Brep: Surface mesh * VDB volume Getting Started --------------- * A simple example of running a structural simulation. Any geometry files used in these examples are available in the ``examples/`` folder. * Import Intact Module :: from PyIntact import * * Create Geometry :: beam_geometry = MeshModel("beam.stl") beam_geometry.instance_id = "simple_cantilever" beam_geometry.refine(0.02) # refine stl for smoother visualization * Create Material (Properties are in MKS unit system) :: material = IsotropicMaterialDescriptor() # steel material.density = 7800.0 # kg/m^3 material.poisson_ratio = 0.3 material.youngs_modulus = 2.1e11 # Pa * Create Restraints and Loads :: # Create a fixed boundary condition at one end of the beam fixed_boundary = FixedBoundaryDescriptor() fixed_boundary.boundary = MeshModel("restraint.stl") # Create a pressure load at the free end of the beam pressure_load = PressureForceDescriptor() pressure_load.units = UnitSystem.MeterKilogramSecond pressure_load.magnitude = 1000 # Applying pressure in Pascals pressure_load.boundary = MeshModel("load.stl") * Setup simulation and solve :: # Setup the simulation scenario descriptor scenario = LinearElasticScenarioDescriptor() scenario.materials = {"Steel": material} scenario.metadata.resolution = 10000 scenario.metadata.units = UnitSystem.MeterKilogramSecond scenario.boundary_conditions = [fixed_boundary, pressure_load] # Associate the material with the model component = MaterialDomain(beam_geometry, "Steel", scenario) assembly = [component] # Initialize and run the simulator simulator = StressSimulator(assembly, scenario) simulator.solve() * Query result and create output file to visualize :: # Query results for stress distribution results = QueryResult(assembly) stress_query = FieldQuery(Field.VonMisesStress) simulator.sample(stress_query, results) # Write results to a file (unit system is saved as a metadata in the vtu file) results.writeVTK("cantilever_beam_results.vtu", UnitSystem.MeterKilogramSecond) How To ====== Creating Simulation Input ------------------------- Units ~~~~~ The default Unit System is **MKS** or ``MeterKilogramSecond``. We allow the following customizations: * Specify custom unit for geometry through the ``Scenario`` unit (default **MKS**) * Specify custom unit for loads (default **MKS**) * Specify custom unit for materials (default **MKS**) +------------------------------+--------------------------+--------------------------+--------------------------+--------------------------+ | Unit System Keywords | Mass | Force | Stress | Length | +==============================+==========================+==========================+==========================+==========================+ | MeterKilogramSecond | kg | N | Pa | m | +------------------------------+--------------------------+--------------------------+--------------------------+--------------------------+ | CentimeterGramSecond | g | dyne | dyne/cm² | cm | +------------------------------+--------------------------+--------------------------+--------------------------+--------------------------+ | MillimeterMegagramSecond | Mg | N | MPa | mm | +------------------------------+--------------------------+--------------------------+--------------------------+--------------------------+ | FootPoundSecond | slug | lbf | lbf/ft² | ft | +------------------------------+--------------------------+--------------------------+--------------------------+--------------------------+ | InchPoundSecond | lbf s²/in | lbf | psi | in | +------------------------------+--------------------------+--------------------------+--------------------------+--------------------------+ Geometry ~~~~~~~~ Geometry is specified through a ``Model`` object * ``MeshModel`` from surface mesh can be defined in two ways. * ``filename`` to define mesh from a file (STL, PLY) * ``instance_id`` a unique name :: mesh_geometry = MeshModel(filename="beam.stl") mesh_geometry.instance_id = "beam" * Constructing mesh face-by-face :: # define an empty mesh mesh_geometry = MeshModel() mesh_geometry.instance_id = "fixed_boundary" # add vertices to the mesh and get the vertex id as output v0 = mesh_geometry.addVertex([0.0, 0.0, 0.0]) v1 = mesh_geometry.addVertex([1.0, 0.0, 0.0]) v2 = mesh_geometry.addVertex([0.0, 1.0, 0.0]) # add faces defined by the vertex id mesh_geometry.addFacet(v0, v1, v2) * Mesh refinement (usually used for finer interpolation of results during visualization). Uses ``MeshModel.refine`` with ``refinement_level`` as an input which is the largest allowed triangle edge as a ratio of the bounding box diagonal. :: # mesh geometry created above is further refined mesh_geometry.refine(refinement_level=0.02) * **VDB Model** * Using a ``.vdb`` file :: vdb_geometry= VDBModel("box.vdb") * In memory using a ``VDB`` object :: import pyopenvdb as vdb #create a vdb double grid object cube = vdb.DoubleGrid() cube.fill(min=(100, 100, 100), max=(199, 199, 199), value=1.0) #create an Intact Model from VDB object vdb_geometry = PyIntact.VDBModel(cube) #create a tesselation of the vdb to get a MeshModel and its list of vertices for sampling mesh = vdb_geometry.tesselate() # mesh is a MeshModel() num_vertices = mesh.vertexCount() print(num_vertices) Material Properties ~~~~~~~~~~~~~~~~~~~ * **Structural** * ``IsotropicMaterialDescriptor`` to create an **Isotropic** material. (Unit is ``MeterKilogramSecond``) * ``density`` is the material density * ``youngs_modulus`` is the elastic modulus * ``poisson_ratio`` is the Poisson ratio :: steel_structural = IsotropicMaterialDescriptor() # Set the density to 7845 kg/m^3, modulus to 200 GPa and poison ratio to 0.29 steel_structural.density = 7845 steel_structural.youngs_modulus = 200e9 steel_structural.poisson_ratio = 0.29 * ``OrthotropicMaterialDescriptor`` to create an **Orthotropic** material. Orthotropic materials are often used to represent composites and wood, where the material properties along one axis are significantly different from the properties along the other axes. (Unit is ``MeterKilogramSecond``) * ``density`` is the material density * ``Ex`` is the elastic modulus in the material’s x-direction * ``Ey`` is the elastic modulus in the material’s y-direction * ``Ez`` is the elastic modulus in the material’s z-direction * ``Gxy`` is the shear modulus in the xy plane * ``Gxz`` is the shear modulus in the xz plane * ``Gyz`` is the shear modulus in the yz plane * ``vxy`` is the poisson ratio in the xy plane * ``vxz`` is the poisson ratio in the xx plane * ``vyz`` is the poisson ratio in the yz plane * ``transform`` is the optional material transform to arbitrarily rotate the axes along which the material properties are defined. Specified as a list of 9 floating point numbers corresponding to a 3x3 matrix. Note the first 3 floating point numbers corresponds to the first row. :: # Create an orthotropic material (Red Pine) material = OrthotropicMaterialDescriptor() material.density = 460.0 # kg/m³ material.Ex = 11.2e9 # Pa material.Ey = 492.8e6 # Pa material.Ez = 985.6e6 # Pa material.Gxy = 907.2e6 # Pa material.Gxz = 1.08e9 # Pa material.Gyz = 123.2e6 # Pa material.vxy = 0.315 material.vxz = 0.347 material.vyz = 0.308 material.transform = [1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0] # optional in this case * **Thermal** * ``ThermalMaterialDescriptor`` to create a thermal material. (Default unit is ``MeterKilogramSecond``) * ``density`` is the material density * ``conductivity`` is the conductivity of the material * ``specific_heat`` is the specific heat of the material :: thermal_material = ThermalMaterialDescriptor() # Set the density, conductivity, and specific heat of the material thermal_material.density = 2700.0 # kg/m^3 thermal_material.conductivity = 170 # W/(m·K) thermal_material.specific_heat = 500 # J/(kg·K) Component and Assembly ~~~~~~~~~~~~~~~~~~~~~~ * A **Component** is created as a ``MaterialDomain`` and takes three inputs * A ``Model`` object (``MeshModel`` or ``VDBModel`` specifically) * A material name * A ``ScenarioDescriptor`` object :: bar1_structural = MaterialDomain(mesh_geometry, "Steel", structural_scenario) bar2_structural = MaterialDomain(vdb_geometry, "Steel", structural_scenario) * **Assembly** is a ``list`` of components that are **bonded** together :: structural_assembly = [bar1_structural, bar2_structural] Boundary Conditions ~~~~~~~~~~~~~~~~~~~ * **Restraints** * **Fixed Boundary** A “Fixed Boundary” restraint fixes the selected geometry in all directions. The fixed boundary only has one input: * the ``boundary`` surface to be fixed :: fixed = FixedVectorDescriptor() # Set the geometry "restraint.ply" to be fixed fixed.boundary = MeshModel("restraint.ply") * **Fixed Vector** A “Fixed Vector” restraint allows for each direction to be optionally set to a specified displacement value, 0 being fixed and ‘None’ being un-restrained. Note that a structural problem must have all three directions restrained somewhere to be valid. Also, note that this boundary condition **doesn’t take custom units** and it’s always in the units of the scenario metadata. A fixed vector requires four inputs: * the ``boundary`` surface to be fixed * the ``x_value`` to set the displacement for the x-axis direction (optional) * the ``y_value`` to set the displacement for the y-axis direction (optional) * the ``z_value`` to set the displacement for the z-axis direction (optional) :: fixed_vector = FixedVectorDescriptor() # Set the geometry "restraint.stl" to be restrained fixed_vector.boundary = MeshModel("restraint.stl") # Set the displacements in each direction fixed_vector.x_value = 0.2 # X-direction displacement of 0.2 m fixed_vector.y_value = None # un-restrained in the Y-direction at this surface fixed_vector.z_value = 0.0 # fixed Z-direction * **Sliding Restraint** A “Sliding Restraint” allows for motion tangential to a specified surface, but fixes motion normal to the specified surface. A sliding restraint only has one input: * the ``boundary`` surface for the sliding restraint condition :: sliding_boundary = SlidingBoundaryDescriptor() # Set the geometry "restraint.stl" for the sliding restraint condition sliding_boundary.boundary = MeshModel("restraint.stl") * **Structural Loads** (``TractionDescriptor``) * **Vector Force** “Vector Force” load is a surface load applied to a face in a specified direction. An example of this load is pressing on the top of a book to push it across a table. A vector load requires four inputs: * the ``boundary`` surfaces where the load is applied * the ``direction`` vector of the force * the ``UnitSystem`` that applies to the magnitude * the ``magnitude`` of the force. :: load = VectorForceDescriptor() # Set the geometry "load.stl" the vector load is applied to load.boundary = MeshModel("load.stl") # Set the vector load direction to be in the -Z direction with a magnitude of 100 lbf load.direction = [0, 0, -1] load.units = UnitSystem.InchPoundSecond load.magnitude = 100 # lbf * **Torque** “Torque” load is a surface load that applies a twisting force around an axis. The direction of the torque is determined using the right-hand rule: using your right hand, point your thumb in the direction of the axis. A positive torque value applies a torque acting in the direction the fingers of your right hand would wrap around the axis. The torque load is applied among the load faces with a distribution that varies linearly from zero at the axis. A Torque load requires five inputs: * the ``boundary`` surfaces where the load is applied * the ``axis`` of rotation about which the torque acts * ``origin`` is the starting point of the axis of rotation * the ``UnitSystem`` that applies to the magnitude * ``magnitude`` of the torque. :: torque_load = TorqueForceDescriptor() # Set the geometry "load.stl" the torque load is applied to torque_load.boundary = MeshModel("load.stl") # Set the axis of the torque axis torque_load.origin = (10, 1, 1) torque_load.axis = [1, 0, 0] # Set the torque to be about the +X (right-hand rule) torque_load.units = UnitSystem.MeterKilogramSecond torque_load.magnitude = 10 # Set the torque magnitude to 10 N*m * **Pressure** A “Pressure” load is a surface load specified in terms of force per unit area. Positive pressures ‘push’ into the surface, and negative pressures ‘pull’. A Pressure Load requires three inputs: * the ``boundary`` surfaces where the load is applied * the ``UnitSystem`` that applies to the magnitude * the ``magnitude`` of the pressure. :: pressure_load = PressureForceDescriptor() # Set the geometry "load.stl" the pressure load is applied to pressure_load.boundary = MeshModel("load.stl") # Set the pressure magnitude to 10 Pa pressure_load.units = UnitSystem.MeterKilogramSecond pressure_load.magnitude = 10 * **Bearing Force** A “Bearing Force” is a surface load applied to a (typically) cylindrical face to approximate the effects of a shaft pressing against the side of a hole. The applied force gets converted to a varying pressure distribution on the portion of the face experiencing compressive pressure. The pressure distribution is computed automatically to achieve the specified overall bearing force. A Bearing Force requires three inputs: * the ``boundary`` surfaces where the load is applied * the ``direction`` vector of the bearing force * the ``UnitSystem`` that applies to the magnitude * the ``magnitude`` of the force :: bearing_load = BearingForceDescriptor() # Set the geometry "load.stl" the bearing load is applied to bearing_load.boundary = MeshModel("load.stl") # Set the loading direction to be in the -Z bearing_load.direction = [0, 0, -1] # Set the magnitude of the load to 100 N bearing_load.units = UnitSystem.MeterKilogramSecond bearing_load.magnitude = 100 * **Thermal Loads** * **Fixed Boundary (Fixed Temperature)** A “Fixed Boundary” load fixes the selected geometry to a specified temperature when the ``value`` is specified and non-zero. The fixed boundary has two inputs: * the ``boundary`` surface to be fixed * the ``value`` of the temperature to fix at the boundary surface :: fixed = FixedBoundaryDescriptor() # Set the geometry "fixed_temp.ply" to set a fixed temperaure of 320 K fixed_boundary.boundary = MeshModel("fixed_temp.ply") fixed_boundary.value = 320 # Kelvin * **Convection** A “Convection” load specifies the transfer of heat from a surrounding medium. A thermal convection boundary condition requires three inputs: * the ``boundary`` surface(s) where the convection is applied. * the heat transfer ``coefficient`` * the ``environment_temperature`` of the surrounding medium * the ``units`` (default MKS) :: # Create an instance of ConvectionDescriptor convection = ConvectionDescriptor() convection.units = UnitSystem.MeterKilogramSecond # Set the geometry "face.ply" the convection is applied to convection.boundary = MeshModel("face.ply") # Set the heat transfer coefficient to 25 W/m^2K and environment temperature convection.coefficient = 25 # W/m^2K convection.environment_temperature = 300 # Kelvin * **Surface Flux** Surface “Thermal or Heat Flux” specifies the heat flow per unit of surface area. A surface thermal flux requires two inputs: * the ``boundary`` surface(s) where the flux is applied. * the ``magnitude`` of the heat flux * the ``units`` (default MKS) :: # Create an instance of ConstantFluxDescriptor flux = ConstantFluxDescriptor() flux.units = UnitSystem.MeterKilogramSecond # Set the geometry "face.ply" the constant flux is applied to flux.boundary = MeshModel("face.ply") # Set the flux magnitude to 500 W/m^2 flux.magnitude = 500 # W/m^2 * **Body Loads/Internal Conditions** Add body load to the scenario as shown below. (see Scenario Setup section for more details) :: scenario.internal_conditions = [rotational_load, gravity_load] * **Body Load (Linear Acceleration, Gravity)** Body loads comprise forces that are distributed over a solid volume. They are specified by entering the components of a linear acceleration vector field in which the body is immersed. The material in the body will tend to be pulled in the direction of the acceleration vector. The inputs to the body load are: * the ``direction`` vector of the acceleration field * the ``UnitSystem`` that applies to the magnitude * the ``magnitude`` of acceleration :: # Example for a "gravity load" body_load = BodyLoadDescriptor() # Set the direction vector to be downward (-Z) body_load.direction = [0, 0, -1] # Set the magnitude to 9.80655 m/s^2 body_load.units = UnitSystem.MeterKilogramSecond body_load.magnitude = 9.80665 * **Rotational Load** Rotational body loads simulate the effect of a body rotating around an axis. Two contributions are considered in a rotational body load: angular velocity and angular acceleration. The angular velocity term simulates the centrifugal effects that tend to throw a body’s material away from the axis of rotation. The angular acceleration term simulates the effect of a rotational acceleration field around the axis of rotation. A positive angular acceleration tends to drag the body’s material in the positive rotational direction according to the right-hand rule. A rotational body load has 4 inputs: * ``origin`` point for the axis of rotation * vector defining the ``axis`` of rotation * ``angular_velocity`` (in radians/sec) * ``angular_acceleration`` (in radians/sec²) :: rotational_load = RotationalLoadDescriptor() # Set the origin at the coordinate system origin rotational_load.origin = (0, 0, 0) # Set the axis of rotation about the y-axis rotational_load.axis = [0, 1, 0] # Set angular velocity to 10 rad/s and angular acceleration to 0.5 rad/s^2 rotational_load.angular_velocity = 10 rotational_load.angular_acceleration = 0.5 * **Thermal Loads** * **Constant Heat** A “Constant Heat” or body heat flux load applies uniform heat generation over a specified volume. Constant heat flux has 2 inputs * the ``instance_id`` of the components which are producing heat flux * the ``magnitude`` of the body heat flux * the ``units`` (default MKS) :: constant_heat = ConstantHeatDescriptor() constant_heat.units = UnitSystem.MeterKilogramSecond # Set the heat generation to -200,000 W/m^3 for a beam component constant_heat.instance_id = "beam" constant_heat.magnitude = -200000.0 # W/m^3 Simulation Scenario Setup and Solution ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ * **Scenario Setup** using ``ScenarioDescriptor`` * **Linear Elasticity** scenario is created using ``LinearElasticScenarioDescriptor``. It takes the following inputs: * ``boundary_conditions``, a list of boundary conditions * ``internal_conditions``, a list of internal conditions (body loads) * **Scenario Metadata,** ``metadata``, a set of solver parameters *to control the accuracy and speed of the simulation*: * ``resolution`` is the target number of finite elements. An iterative process determines a ``cell_size`` that achieves approximately the specified number of elements. You can directly specify the ``cell_size`` instead of ``resolution``. Note that decreasing ``cell_size`` can quickly result in large numbers of finite elements and long solve times. ``resolution`` is recommended for most cases. * ``units`` sets the unit of the scenario (geometry and results will be in this unit) **for example, a geometry in ‘MKS’ that is 6 m long would be 6 mm long when set to ‘MMS’** * ``basis_order`` is the type of finite element used. The default is 1, which would be linear elements. 2 is for quadratic elements. * ``solver_type`` is the solver type used in the simulation with the following options: * ``MKL_PardisoLDLT`` (**default**) is the direct solver and is typically faster for lower resolution (< 200K cells on 32 GB memory), but gets slower and consumes more memory at higher resolutions * ``AMGCL_amg_rigid_body`` is the iterative solver which is typically faster at high resolution. :: scenario = LinearElasticScenarioDescriptor() scenario.materials = {"Aluminum 6061-T6": material} scenario.boundary_conditions = [fixed_boundary, vector_load, torque_load] scenario.internal_conditions = [body_load] scenario.metadata.resolution = 10000 scenario.metadata.units = UnitSystem.MeterKilogramSecond #Optional settings scenario.metadata.basis_order = 2 #2 = quadratic elements scenario.metadata.solver_override = Solver.AMGCL_amg_rigid_body #iterative solver * Similarly, **Modal** scenario is created using ``ModalScenarioDescriptor`` which takes ``boundary_conditions``, a list of boundary conditions, as input. * The ``metadata`` is the same as for the linear elastic scenario except there is necessary input ``metadata.desired_eigenvalues`` that takes in an integer for the number of eigenvalues. :: # Setup the modal simulation scenario modal_scenario = ModalScenarioDescriptor() modal_scenario.materials = {"Aluminum 6061-T6": material} modal_scenario.boundary_conditions = [fixed_boundary] modal_scenario.metadata.resolution = 10000 modal_scenario.modal_metadata.desired_eigenvalues = 10 #Optional settings modal_scenario.metadata.basis_order = 2 #2 = quadratic elements modal_scenario.metadata.solver_override = Solver.AMGCL_amg_rigid_body #iterative solver * Similarly, a thermal scenario is created using ``StaticThermalScenarioDescriptor`` which takes a list of ``boundary_conditions`` and a list of ``internal_conditions`` as input. Note that the ``material`` must be a thermal material for this scenario type. * The ``metadata`` is the same as for the linear elastic scenario with additional input ``metadata.environment_temperature`` which specifies the temperature of the environment. :: # Setup the thermal simulation scenario thermal_scenario = StaticThermalScenarioDescriptor() thermal_scenario.materials = {"Aluminum 6061-T6": material} thermal_scenario.boundary_conditions = [fixed_temp, convection] thermal_scenario.internal_conditions = [constant_heat] thermal_scenario.thermal_metadata.environment_temperature = 0.0 thermal_scenario.metadata.resolution = 10000 #Optional settings thermal_scenario.metadata.basis_order = 2 #2 = quadratic elements thermal_scenario.metadata.solver_override = Solver.AMGCL_amg #iterative solver * **Finalize Simulation Setup and Execute** * **Stress Simulation** is created using ``StressSimulator``. It takes the following inputs: * An ``assembly`` or the list of ``MaterialDomains`` * ``LinearElasticScenarioDescriptor`` that describes the simulation scenario :: # Initialize and run the linear elastic scenario simulator = StressSimulator(assembly, scenario) simulator.solve() * **Modal Simulation** is created using ``ModalSimulator``. It takes the following inputs: * An ``assembly`` or the list of ``MaterialDomains`` * ``ModalScenarioDescriptor`` that describes the simulation scenario :: # Intialize and run the modal scenario modal_simulator = ModalSimulator(assembly, modal_scenario) modal_simulator.solve() Query and Result Output ----------------------- * Create a ``QueryResult`` object that specifies the domain to be sampled. It takes a single ``component`` or an ``assembly`` as input :: assembly = [component_1, component_2] results = QueryResult(component_1) # sample only the defined component results = QueryResult(assembly) # sample on the full assembly * Specify a **query class & type** that defines what quantity you are querying. The query classes and types available are as follows : * **Global Query** is for querying quantities defined for the entire domain/structure. The two ``GlobalQueryType`` subclasses available are ``Frequency`` & ``Compliance`` :: # GlobalQuery example # Create a GlobalQuery to query the first mode results = QueryResult(assembly) frequency_query = GlobalQuery(GlobalQueryType.Frequency, DiscreteIndex(0)) # Sample the simulation for the specified query freq_1 = simulator.sample(mode1_query, results) print(freq_1.get(0,0)) # Print out the first mode in Hz # Example use to get the first 10 modes # Note, the simulator would need this: # modal_scenario.modal_metadata.desired_eigenvalues = 10 for i in range(10): # Create the query for the current mode query = GlobalQuery(GlobalQueryType.Frequency, DiscreteIndex(i)) # Sample the query using the simulator r = simulator.sample(query, results) # Append the value obtained from result.get(0, 0) to a list values.append(r.get(0, 0)) print(values) # list of first 10 modes in Hz * **Field Query** is for querying field quantities defined at any point in the domain. This depends on the physics type. * For **Stress Simulation**, the following ``FieldType`` inputs are available * ``Displacement`` tuple of dimension 3 for each point :: [x-displacement, y-displacement, z-displacement] * ``Strain`` and ``Stress`` are each a tuple of dimension 6 for each point :: [stress_xx, stress_yy, stress_zz, stress_yz, stress_xz, stress_xy] [strain_xx, strain_yy, strain_zz, strain_yz, strain_xz, strain_xy] * ``TopologicalSensitivity``, ``StrainEnergyDensity``, and ``VonMisesStress`` are each a tuple of dimension 1 for each point * **Field Query Interface usage** :: # FieldQuery example # Create a FieldQuery to query displacement results = QueryResult(assembly) displacement_query = FieldQuery(f=Field.Displacement) # Create a FieldQuery with optional input to query the y-component y_displacement_query = FieldQuery(f=Field.Displacement, component=1) # FieldQuery with optional input to query the norm of displacement displacement_magnitude_query = FieldQuery(f=Field.Displacement, norm=True) # Sample the simulation with the given query to create a VectorArray # * more info on VectorArrays are provided in the subsequent section * displacement_VectorArray1 = simulator.sample(displacement_query, results) displacement_VectorArray2 = simulator.sample(y_displacement_query, results) displacement_VectorArray3 = simulator.sample(displacement_magnitude_query, results) displacement1 = displacement_VectorArray1.get(i, j); # jth component of the displacement (0 = X, 1 = Y, or 2 = Z) at the ith sampling point displacement2 = displacement_VectorArray2.get(1, 0); # y-displacement value at the second sampling point displacement3 = displacement_VectorArray3.get(1, 0); # displacement magnitude at the second sampling point * For **Modal Simulation**, the only ``FieldType`` available is ``Displacement``, which corresponds to the mode shape information and requires an index for the specified mode. :: # Create a field query for the displacement/mode shape of mode 1 mode_num = 0 # first mode modal_field_query = FieldQuery(f=Field.Displacement, scheme=DiscreteIndex(mode_num)) Sample and Export Results ~~~~~~~~~~~~~~~~~~~~~~~~~ * **Sampling results in-memory** * Methods to sample results of the simulation and store/print ``VectorArray`` results from the queries created above. :: # Query results results = QueryResult(assembly) stress_query = FieldQuery(f=Field.Stress) # a field is needed to create a field query stress_VectorArray = simulator.sample(stress_query, results) * This sampling is stored in a ``VectorArray`` of size ``n_tuple`` by ``dimension``, ``n_tuple`` is the number of points sampled in the component/assembly and dimension depends on the query type. For example, ``FieldQuery(Field.Displacement)`` has a dimension of 3 for each displacement component, thus ``VectorArray.get(i, 2)`` would get the z-displacement of point ``i``. :: # Get the number of tuples (vectors) in the VectorArray (one per point results are sampled at) n_tuples = stress_VectorArray.n_tuples() # Get the dimension of each vector - in this case 6, one value for each stress component dimension = stress_VectorArray.dimension() stress_xx = stress_VectorArray.get(i, 0); # stress_xx at i-th sample point stress_yy = stress_VectorArray.get(i, 1); # stress_yy at i-th sample point stress_ij = stress_VectorArray.get(i, j); # jth stress component (0, 1, 2, 3, 4, or 5) of the ith tuple * **Export results to a file**: Once the query is created, ``results`` can be used to write a ``.vtu`` file which contains **all** solution fields via ``results.writeVTK`` which has two inputs: * file name ``"*.vtu"`` * ``UnitSystem`` is an attribute in the ``*.vtu`` ( Note that this argument is just metadata in the ``*.vtu`` file. The result magnitudes are in the unit of the scenario regardless of the unit system specified. ) :: # Write results to a file (unit system is saved as a metadata in the vtu file) results = QueryResult(assembly) stress_query = FieldQuery(Field.VonMisesStress) simulator.sample(stress_query, results) # Write results to a file (unit system is saved as a metadata in the vtu file) results.writeVTK("cantilever_beam_results.vtu", UnitSystem.MeterKilogramSecond) Sampling on Custom Geometries ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ * Another important functionality is the ability to specify a point set to sample on. This can be done by creating a new ``sampling_component`` consisting of vertices at the desired sampling locations. **Note** that the ``sampling_component`` is not to be included in the assembly that is being simulated. :: # Create an empty MeshModel to store vertices for desired sampling locations sampling_geometry = MeshModel() # Add the desired sampling point locations and store their ID for later use test_index1 = sampling_geometry.addVertex([6.0, 0.0, 0.0]) test_index2 = sampling_geometry.addVertex([3.0, 0.0, 0.0]) sampling_component = MaterialDomain(sample_geometry, "material_name", scenario) # Define the assembly to be simulated (note that the sampling_component is not included) component = MaterialDomain(component_geometry, "material_name", scenario) assembly = [component] simulator = StressSimulator(assembly, scenario) simulator.solve() displacement_query = FieldQuery(f=Field.Displacement) sampled_results = QueryResult(sampling_component) # sample only the custom point set sampled_VectorArray = simulator.sample(displacement_query, sampled_results) displacement_index1 = sampled_VectorArray.get(test_index1, 0) # x-displacement at above specified index1 displacement_index2 = sampled_VectorArray.get(test_index2, 1) # y-displacement at above specified index2 Installation ============ Prerequisites ------------- PyIntact can be installed when running Python 3.10 or 3.11, and requires the ``openmpi-devel`` package to be installed. Installation ------------ Install PyIntact, the python interface to Intact.Simulation, by downloading the Python wheel and then using ``pip`` to install, as in the following example: :: pip install PyIntact.whl PyIntact API ============ **class PyIntact.Geometry(self: PyIntact.Geometry, value: int)** Members: VDB Mesh ``property name`` **class PyIntact.Model** ``property instance_id`` The instance ID of the model. We enforce that this must be unique per run. **mass(self: PyIntact.Model) -> float** Returns the mass of the model ``property model_name`` This is used for debugging purposes. It need not be filled out in general. **projectField(self: PyIntact.Model, field: Callable[[Annotated[list[float], FixedSize(3)]], Intact::Tensor]) -> None** Sets tensors according to a user-provided callback Parameter ``field``: The field **setDensities(self: PyIntact.Model, densities: list[float]) -> None** Sets tensors as 1x1-dimensional data, AKA densities. This is a convenience wrapper around *setTensors*. Parameter ``densities``: The densities **setTensors(self: PyIntact.Model, tensors: list[Intact::Tensor]) -> None** Set the tensors for all voxels. Parameter ``tensors``: The tensors for all voxels. **tensor(self: PyIntact.Model, id: int) -> Intact::Tensor** Returns the tensor at the given index Parameter ``id``: The id (typically point ID) for the tensor **class PyIntact.MeshModel(*args, **kwargs)** Bases: `Model <#PyIntact.Model>`_ Overloaded function. 1. __init__(self: PyIntact.MeshModel) -> None Construct an empty MeshModel. 1. __init__(self: PyIntact.MeshModel, other: PyIntact.MeshModel) -> None Construct a MeshModel from an existing MeshModel. Parameter ``other``: A MeshModel. 1. __init__(self: PyIntact.MeshModel, filename: str) -> None Load MeshModel from file. Parameter ``filename``: File must be a VTU, PLY, or STL file. **addFacet(self: PyIntact.MeshModel, vertex_1: int, vertex_2: int, vertex_3: int) -> int** Add a new facet to the MeshModel. Parameter ``vertex_1``: Index of the first vertex in the facet. Parameter ``vertex_2``: Index of the second vertex in the facet. Parameter ``vertex_3``: Index of the third vertex in the facet. :Returns: The index of the new facet. **addVertex(self: PyIntact.MeshModel, vertex: Annotated[list[float], FixedSize(3)]) -> int** Add new vertex to the MeshModel. Parameter ``vertex``: The position of the vertex. :Returns: The index of the new vertex for use with facets. **facet(self: PyIntact.MeshModel, facet_index: int) -> Annotated[list[int], FixedSize(3)]** Get a specific facet. Parameter ``facet_index``: The facet identifier. :Returns: The facet. **facetCount(self: PyIntact.MeshModel) -> int** Get the number of facets in the MeshModel. **refine(*args, **kwargs)** Overloaded function. 1. refine(self: PyIntact.MeshModel, threshold: float) -> None Refine the input mesh, by splitting every facet edge that is larger than (threshold) * (bounding box size of mesh). Refining the input mesh increases sampling density. Parameter ``threshold``: The fraction of the bounding box size of the mesh. 1. refine(self: PyIntact.MeshModel) -> None Refine the input mesh, by splitting every facet edge that is larger than 2% of bounding box size of mesh. **vertex(self: PyIntact.MeshModel, vertex_index: int) -> Annotated[list[float], FixedSize(3)]** Get a specific vertex. Parameter ``vertex_index``: The vertex identifier. :Returns: The vertex. **vertexCount(self: PyIntact.MeshModel) -> int** Get the number of vertices in the MeshModel. **writePLY(self: PyIntact.MeshModel, output_filename: str) -> None** Write the MeshModel to file in PLY format. Parameter ``output_filename``: The output filename. **writeVTU(self: PyIntact.MeshModel, output_filename: str) -> None** Write the MeshModel to file in VTU (VTK unstructured grid) format. Parameter ``output_filename``: The output filename. **class PyIntact.Tensor(self: PyIntact.Tensor, arg0: list[float], arg1: int, arg2: int)** **cols(self: PyIntact.Tensor) -> int** **raw_data(self: PyIntact.Tensor) -> list[float]** **rows(self: PyIntact.Tensor) -> int** **class PyIntact.VDBModel(*args, **kwargs)** Bases: `Model <#PyIntact.Model>`_ Represents a model based on an openVDB volume. Overloaded function. 1. __init__(self: PyIntact.VDBModel, filename: str) -> None Load a VDBModel from a file. Parameter ``filename``: The VDB file. 1. __init__(self: PyIntact.VDBModel, grid: openvdb::v7_2::Grid, 4u>, 5u> > > >) -> None Construct a VDBModel an existing DoubleGrid. Parameter ``grid``: The existing DoubleGrid. **tessellate(self: PyIntact.VDBModel) -> `PyIntact.MeshModel <#PyIntact.MeshModel>`_** Create a mesh from the VDBModel, using the OpenVDB volumeToMesh method. :Returns: A ``MeshModel``. **class PyIntact.MaterialDomain(self: PyIntact.MaterialDomain, model: PyIntact.Model, material_name: str, scenario_descriptor: Intact::AbstractScenarioDescriptor)** Construct a MaterialDomain from a model, material name, and scenario descriptor. This constructor ensures the material units are consistent with scenario units and should be used. Parameter ``model``: The model Parameter ``material_name``: The material name Parameter ``scenario_descriptor``: The scenario descriptor **class PyIntact.Field(self: PyIntact.Field, value: int)** Members: Displacement Strain Stress TopologicalSensitivity StrainEnergyDensity VonMisesStress Temperature HeatFlux ``property name`` **class PyIntact.Query** **class PyIntact.FieldQuery(*args, **kwargs)** Bases: `Query <#PyIntact.Query>`_ This class describes a query over a field. This query can be for a specific field component, or the norm of the field. Overloaded function. 1. __init__(self: PyIntact.FieldQuery, *, f: PyIntact.Field, norm: bool = False) -> None Parameter ``f``: The field to sample Parameter ``norm``: Whether to compute the norm of the field. 1. __init__(self: PyIntact.FieldQuery, *, f: PyIntact.Field, scheme: Intact::IndexScheme, norm: bool = False) -> None Parameter ``f``: The field to sample Parameter ``scheme``: The indexing scheme for the field Parameter ``norm``: Whether to compute the norm of the field. 1. __init__(self: PyIntact.FieldQuery, *, f: PyIntact.Field, component: int) -> None Parameter ``f``: The field to sample Parameter ``component``: The field component to sample 1. __init__(self: PyIntact.FieldQuery, *, f: PyIntact.Field, component: int, scheme: Intact::IndexScheme) -> None Parameter ``f``: The field to sample Parameter ``component``: The field component to sample Parameter ``scheme``: The indexing scheme for the field **class PyIntact.StatisticalQuery(*args, **kwargs)** Bases: `Query <#PyIntact.Query>`_ This class describes a query for a statistic over a field. This query can be for the minimum, maximum or mean value of the field. Overloaded function. 1. __init__(self: PyIntact.StatisticalQuery, *, f: PyIntact.Field, s: PyIntact.StatisticalQueryType, norm: bool = False) -> None Parameter ``f``: The field to sample. Parameter ``s``: The statistic to query for. Parameter ``norm``: Whether to compute the norm of the field. 1. __init__(self: PyIntact.StatisticalQuery, *, f: PyIntact.Field, component: int, s: PyIntact.StatisticalQueryType) -> None Parameter ``f``: The field to sample. Parameter ``component``: The field component to sample. Parameter ``s``: The statistic to query for. 1. __init__(self: PyIntact.StatisticalQuery, *, f: PyIntact.Field, s: PyIntact.StatisticalQueryType, scheme: Intact::IndexScheme, norm: bool) -> None Parameter ``f``: The field to sample. Parameter ``s``: The statistic to query for. Parameter ``scheme``: The indexing scheme for the field Parameter ``norm``: Whether to compute the norm of the field. 1. __init__(self: PyIntact.StatisticalQuery, *, f: PyIntact.Field, component: int, s: PyIntact.StatisticalQueryType, scheme: Intact::IndexScheme) -> None Parameter ``f``: The field to sample. Parameter ``component``: The field component to sample. Parameter ``s``: The statistic to query for. Parameter ``scheme``: The indexing scheme for the field **class PyIntact.IndexScheme** **class PyIntact.DiscreteIndex(self: PyIntact.DiscreteIndex, arg0: int)** Bases: `IndexScheme <#PyIntact.IndexScheme>`_ **class PyIntact.NullIndex** Bases: `IndexScheme <#PyIntact.IndexScheme>`_ **class PyIntact.VectorArray(*args, **kwargs)** This class represents an array of constant-sized vector data. It wraps a flat std::vector for cache performance, but structures the queries so that indexing errors are unlikely. Overloaded function. 1. __init__(self: PyIntact.VectorArray) -> None 2. __init__(self: PyIntact.VectorArray, arg0: int, arg1: int) -> None 3. __init__(self: PyIntact.VectorArray, arg0: list[float], arg1: int) -> None **dimension(self: PyIntact.VectorArray) -> int** **get(self: PyIntact.VectorArray, arg0: int, arg1: int) -> float** Get the component at the given index. Parameter ``index``: The index Parameter ``component``: The component :Returns: Returns the component at the given index value. **n_tuples(self: PyIntact.VectorArray) -> int** Gives the number of tuples stored in this vector array. :Returns: The full size of the raw data divided by the vector m_dimension. **raw_data(self: PyIntact.VectorArray) -> list[float]** **class PyIntact.QueryResult(*args, **kwargs)** This class describes a query result. The buffer contains a list of scalar values, however the amount of data depends on the query type that produced it. Overloaded function. 1. __init__(self: PyIntact.QueryResult, model: PyIntact.MaterialDomain) -> None Constructs a new query result over the given material domain. Note that if the underlying model is a VoxelModel, the data will be sampled at the centroids of each voxel, and the data will be persisted in the same order as the voxels themselves. Otherwise, the data will be sampled at the points of the model and persisted in the same order as the points. 1. __init__(self: PyIntact.QueryResult, assembly: list[PyIntact.MaterialDomain]) -> None Constructs a new query result over the given assembly. Note that if the underlying model is a VoxelModel, the data will be sampled at the centroids of each voxel, and the data will be persisted in the same order as the voxels themselves. Otherwise, the data will be sampled at the points of the model and persisted in the same order as the points. **addData(self: PyIntact.QueryResult, arg0: PyIntact.Query, arg1: PyIntact.VectorArray) -> None** Adds the result of computing *query* to this object. Parameter ``query``: The query being computed. Parameter ``data``: The data result of the query. **getData(self: PyIntact.QueryResult, arg0: PyIntact.Query) -> `PyIntact.VectorArray <#PyIntact.VectorArray>`_** Returns the vector-valued data indexed by the given query. Note that this method can throw an error if the query data has not been precomputed. Use *hasData* to ensure that such data has been computed. Parameter ``query``: The query whose results are to be returned :Returns: The vector-valued output of the query. **hasData(self: PyIntact.QueryResult, arg0: PyIntact.Query) -> bool** Determines if data exists for a given query. Parameter ``query``: The query whose results we are checking on. :Returns: True if data has been computed, False otherwise. **writeVTK(self: PyIntact.QueryResult, filename: str, unit_system: Intact::UnitSystemType) -> None** Writes a vtk file containing the results that have been stored on this object. Parameter ``filename``: The filename to write to. If it does not end in “.vtu”, the “.vtu” will be appended to the filename. Parameter ``unit_system``: The Unitsystem to use when attaching units to the sampled data. **class PyIntact.UnitSystem(self: PyIntact.UnitSystem, value: int)** Members: MeterKilogramSecond CentimeterGramSecond MillimeterMegagramSecond FootPoundSecond InchPoundSecond ``property name`` **class PyIntact.Metadata** ``property basis_order`` The order of basis function approximation to use. We support linear (= 1) and quadratic (= 2). ``property cell_size`` The exact simulation cell size to use. Warning: Setting the cell_size overrides the resolution. ``property integrator_override`` For some physics types, there may be multiple integrator implementations available. This field allows you to override the default solver. ``property max_iterations_multiplier`` max_iterations_multiplier is used to control the maximum number of iterations of a solver from the AMG family of solvers. Use a value less than 1 to reduce the maximum number of iterations and reduce the solver runtime. ``property moment_order`` moment_order is the order of moments computed for integration. Use lower for faster runtime, at an accuracy cost. This number must be at least 2 and can be as high as 6-7. The default is 3. ``property resolution`` resolution is the number of finite elements to compute over. Note that Intact will only approximate the target number of cells. Either this member or the *cell_size* must be set, or else Intact will throw an error. ``property solver_override`` For some integrator types, there are multiple sparse linear solvers available. Please see the user’s manual for information about the available options. The default is AMGCL_amg. ``property tolerance`` The tolerance of the linear solver. Use a smaller tolerance for more accurate simulations and longer runtime. Use a larger tolerance for less accurate simulations and shorter runtime. Default is 1e-8. **class PyIntact.ThermalMetadata** ``property environment_temperature`` The temperature of the environment. **class PyIntact.ModalMetadata** ``property desired_eigenvalues`` The number of eigenvalues or natural frequencies to calculate. **class PyIntact.Solver(self: PyIntact.Solver, value: int)** Members: AMGCL_cg : Iterative conjugate gradient solver without preconditioner. AMGCL_ilu0_cg : Iterative conjugate gradient solver with ilu0 preconditioner. AMGCL_damped_jacobi_cg : Iterative conjugate gradient solver with damped jacobi preconditioner. AMGCL_amg : Iterative AMG based solver with smoothed aggregation. This is the default solver for thermal scenarios. AMGCL_amg_rigid_body : Iterative AMG based solver with rigid body modes based coarsening. This is the default solver for linear elasticity. This solver cannot be used with thermal scenarios. Eigen_SparseLU : Direct solver that uses LU decomposition. MKL_PardisoLDLT : PardisoLDLT is proprietary direct solver provided by Intel through MKL. Use Intact::mklPresent() to check whether this feature is available. MKL_PardisoLLT : PardisoLLT is proprietary direct solver provided by Intel through MKL. Use Intact::mklPresent() to check whether this feature is available. MKL_PardisoLU : PardisoLU is proprietary direct solver provided by Intel through MKL. Use Intact::mklPresent() to check whether this feature is available. ``property name`` **class PyIntact.Integrator(self: PyIntact.Integrator, value: int)** Members: IntactModal IntactLinearElasticity MPCLinearElasticity MPCStaticThermal LinearBuckling ``property name`` **class PyIntact.AbstractBoundaryConditionDescriptor** ``property boundary`` The boundary where the boundary condition is applied. ``property units`` The unit system. If not set, then the default value is ``MeterKilogramSecond``. **class PyIntact.TractionDescriptor** Bases: `AbstractBoundaryConditionDescriptor <#PyIntact.AbstractBoundaryConditionDescriptor>`_ **class PyIntact.VectorForceDescriptor(self: PyIntact.VectorForceDescriptor)** Bases: `TractionDescriptor <#PyIntact.TractionDescriptor>`_ Describes a vector force applied uniformly over the specified boundary. ``property direction`` The direction of the force. ``property magnitude`` The magnitude of the force. **class PyIntact.HydrostaticForceDescriptor(self: PyIntact.HydrostaticForceDescriptor)** Bases: `TractionDescriptor <#PyIntact.TractionDescriptor>`_ Simulates a body submerged in liquid. ``property density`` The density of the liquid in which the body is submerged. ``property height`` The height the liquid rises above the ``z=0`` plane. **class PyIntact.PressureForceDescriptor(self: PyIntact.PressureForceDescriptor)** Bases: `TractionDescriptor <#PyIntact.TractionDescriptor>`_ Describes a pressure load (force per unit area) applied normal to the boundary. ``property magnitude`` The magnitude of the pressure load. **class PyIntact.TorqueForceDescriptor(self: PyIntact.TorqueForceDescriptor)** Bases: `TractionDescriptor <#PyIntact.TractionDescriptor>`_ Describes a torque load applied to the boundary. ``property axis`` The axis of rotation of the torque load. ``property magnitude`` The magnitude of the torque load. ``property origin`` The origin of the axis of rotation of the torque. **class PyIntact.BearingForceDescriptor(self: PyIntact.BearingForceDescriptor)** Bases: `TractionDescriptor <#PyIntact.TractionDescriptor>`_ Describes a load at the contact between two curved surfaces. The contact pressure is not uniform. It is scaled by the angle of the contact relative to the angle of the force. ``property direction`` The direction of the force. ``property magnitude`` The magnitude is the total force applied over the entire surface. **class PyIntact.FixedBoundaryDescriptor(self: PyIntact.FixedBoundaryDescriptor)** Bases: `AbstractBoundaryConditionDescriptor <#PyIntact.AbstractBoundaryConditionDescriptor>`_ This class describes a fixed solution quantity (also known as a Dirichlet boundary condition). In order to set the value to fix the solution, set the ``value`` member. Note that not all scenarios support inhomogeneous Dirichlet conditions. ``property value`` The value of the boundary condition. **class PyIntact.FixedVectorDescriptor(self: PyIntact.FixedVectorDescriptor)** Bases: `AbstractBoundaryConditionDescriptor <#PyIntact.AbstractBoundaryConditionDescriptor>`_ Describes a fixed vector quantity. A fixed vector boundary condition can represent fixing displacement at some non-zero value, or it can represent partial restraints, also known as axis-aligned sliding restraints. Displacement is a three dimensional field, and each displacement component can optionally be specified. Any axis not specified will not be restrained. ``property x_value`` The value of the boundary condition in the x-direction. ``property y_value`` The value of the boundary condition in the y-direction. ``property z_value`` The value of the boundary condition in the z-direction. **class PyIntact.SlidingBoundaryDescriptor(self: PyIntact.SlidingBoundaryDescriptor)** Bases: `AbstractBoundaryConditionDescriptor <#PyIntact.AbstractBoundaryConditionDescriptor>`_ Describes a sliding boundary condition, which prevents displacement in the direction normal to the specified surface, and permits displacement tangential to the surface. **class PyIntact.ConvectionDescriptor(self: PyIntact.ConvectionDescriptor)** Bases: `AbstractBoundaryConditionDescriptor <#PyIntact.AbstractBoundaryConditionDescriptor>`_ Describes the convective transfer of heat from a surrounding medium. ``property coefficient`` The coefficient of convection. ``property environment_temperature`` The temperature of the surrounding medium. **class PyIntact.ConstantFluxDescriptor(self: PyIntact.ConstantFluxDescriptor)** Bases: `AbstractBoundaryConditionDescriptor <#PyIntact.AbstractBoundaryConditionDescriptor>`_ Describes a constant flux per unit surface area. ``property magnitude`` The magnitude of the flux. **class PyIntact.AbstractInternalConditionDescriptor** ``property instance_id`` The instance id of the ``Model`` that this internal condition is applied to. If not set, then this internal condition is applied to every instance. ``property units`` The unit system. If not set, then the default value is ``MeterKilogramSecond``. **class PyIntact.BodyLoadDescriptor(self: PyIntact.BodyLoadDescriptor)** Bases: `AbstractInternalConditionDescriptor <#PyIntact.AbstractInternalConditionDescriptor>`_ Describes body load caused by linear acceleration. ``property direction`` The direction of the acceleration. ``property magnitude`` The magnitude of the acceleration. **class PyIntact.RotationalLoadDescriptor(self: PyIntact.RotationalLoadDescriptor)** Bases: `AbstractInternalConditionDescriptor <#PyIntact.AbstractInternalConditionDescriptor>`_ Describes a rotational load on a body. ``property angular_acceleration`` The angular acceleration. ``property angular_velocity`` The angular velocity. ``property axis`` The axis of rotation of the rotational load. ``property origin`` The origin of the axis of rotation of the rotational load. **class PyIntact.ConstantHeatDescriptor(self: PyIntact.ConstantHeatDescriptor)** Bases: `AbstractInternalConditionDescriptor <#PyIntact.AbstractInternalConditionDescriptor>`_ Describes the uniform generation of heat within a body. ``property magnitude`` The magnitude of the body heat flux. **class PyIntact.AbstractMaterialDescriptor** An abstract description of a material. ``property density`` The material density. ``property units`` The unit system. If not set, then the default value is ``MeterKilogramSecond``. **class PyIntact.IsotropicMaterialDescriptor(self: PyIntact.IsotropicMaterialDescriptor)** Bases: `AbstractMaterialDescriptor <#PyIntact.AbstractMaterialDescriptor>`_ Describes an isotropic material, which has the same behavior no matter the direction of the forces applied. ``property compressive_strength`` The ultimate strength in compression of the material. ``property poisson_ratio`` The Poisson ratio of the material. ``property tensile_strength`` The ultimate strength in tension of the material. ``property yield_strength`` The yield strength of the material. ``property youngs_modulus`` The Young’s modulus of the material. **class PyIntact.OrthotropicMaterialDescriptor(self: PyIntact.OrthotropicMaterialDescriptor)** Bases: `AbstractMaterialDescriptor <#PyIntact.AbstractMaterialDescriptor>`_ ``property Ex`` The elastic modulus in the x-direction ``property Ey`` The elastic modulus in the y-direction ``property Ez`` The elastic modulus in the z-direction ``property Gxy`` The shear modulus in the xy plane. ``property Gxz`` The shear modulus in the xz plane. ``property Gyz`` The shear modulus in the yz plane. ``property transform`` A transformation matrix to change the material orientation. ``property vxy`` The poisson ratio in the xy plane. ``property vxz`` The poisson ratio in the xz plane. ``property vyz`` The poisson ratio in the yz plane. **class PyIntact.ThermalMaterialDescriptor(self: PyIntact.ThermalMaterialDescriptor)** Bases: `AbstractMaterialDescriptor <#PyIntact.AbstractMaterialDescriptor>`_ Describes the thermal properties of an isotropic material. ``property conductivity`` The thermal conductivity. ``property expansion_coefficient`` The coefficient of thermal expansion (optional). ``property specific_heat`` The specific heat of the material. **class PyIntact.AbstractScenarioDescriptor** ``property boundary_conditions`` The collection of boundary conditions to be applied. ``property materials`` The mapping of a material identifier, usually the material name, to the material. ``property metadata`` Metadata about the simulation to be performed. **class PyIntact.LinearElasticScenarioDescriptor(self: PyIntact.LinearElasticScenarioDescriptor)** Bases: `AbstractScenarioDescriptor <#PyIntact.AbstractScenarioDescriptor>`_ Describes a linear elastic scenario to be simulated. ``property internal_conditions`` The collection of the internal conditions, such as a gravity load, to be applied. **class PyIntact.ModalScenarioDescriptor(self: PyIntact.ModalScenarioDescriptor)** Bases: `AbstractScenarioDescriptor <#PyIntact.AbstractScenarioDescriptor>`_ **class PyIntact.StaticThermalScenarioDescriptor(self: PyIntact.StaticThermalScenarioDescriptor)** Bases: `AbstractScenarioDescriptor <#PyIntact.AbstractScenarioDescriptor>`_ **class PyIntact.Simulator** Generic Simulation Manager **cellComplex(self: PyIntact.Simulator) -> list[tuple[Annotated[list[float], FixedSize(3)], Annotated[list[float], FixedSize(3)]]]** Get the simulation grid. :Returns: A vector of pairs of points that represent the minimum and maximum of the bounding box of each cell in the simulation grid. **cellCount(self: PyIntact.Simulator) -> int** Get the number of cells in the simulation grid. **cellSize(self: PyIntact.Simulator) -> float** Get the size of each cell in the simulation grid. **sample(*args, **kwargs)** Overloaded function. 1. sample(self: PyIntact.Simulator, query: PyIntact.Query, result: PyIntact.QueryResult) -> PyIntact.VectorArray Sample the given query over the domain provided. Parameter ``query``: The result to query for. Parameter ``result``: The query result. :Returns: The values of the query. 1. sample(self: PyIntact.Simulator, sample_points: list[Annotated[list[float], FixedSize(3)]], material_name: str, field_queries: list[PyIntact.FieldQuery]) -> list[PyIntact.VectorArray] Sample at the given points, producing results for each query Parameter ``sample_points``: The sample points Parameter ``material_name``: Name of the material associated with points Parameter ``field_queries``: The field queries to be answered :Returns: The values of the field queries at the query points, in order **solve(self: PyIntact.Simulator) -> PyIntact.SolutionStats** Solve the simulation that has been constructed. :Returns: The number of iterations performed by the solver and the error. **class PyIntact.StressSimulator(*args, **kwargs)** Bases: `Simulator <#PyIntact.Simulator>`_ Simulation Manager for Stress Scenarios Overloaded function. 1. __init__(self: PyIntact.StressSimulator, assembly: list[PyIntact.MaterialDomain], descriptor: PyIntact.LinearElasticScenarioDescriptor) -> None Construct a linear elastic scenario Parameter ``assembly``: A collection of MaterialDomains. Parameter ``descriptor``: A LinearElasticScenarioDescriptor that describes the simulation. 1. __init__(self: PyIntact.StressSimulator, arg0: list[PyIntact.MaterialDomain], arg1: PyIntact.LinearElasticScenarioDescriptor, arg2: Intact::StaticThermalSimulator) -> None Construct a linear elastic scenario Parameter ``assembly``: A collection of MaterialDomains. Parameter ``descriptor``: A LinearElasticScenarioDescriptor that describes the simulation. **class PyIntact.ModalSimulator(self: PyIntact.ModalSimulator, arg0: list[PyIntact.MaterialDomain], arg1: PyIntact.ModalScenarioDescriptor)** Bases: `Simulator <#PyIntact.Simulator>`_ Simulation Manager for Modal Scenarios Construct a modal elastic scenario. Parameter ``assembly``: A collection of MaterialDomains. Parameter ``descriptor``: A ModalScenarioDescriptor that describes the simulation. **class PyIntact.StaticThermalSimulator(self: PyIntact.StaticThermalSimulator, arg0: list[PyIntact.MaterialDomain], arg1: PyIntact.StaticThermalScenarioDescriptor)** Bases: `Simulator <#PyIntact.Simulator>`_ Simulation Manager for Static Thermal Scenarios Construct a static thermal scenario Parameter ``assembly``: A collection of MaterialDomains. Parameter ``descriptor``: A StaticThermalScenarioDescriptor that describes the simulation. **PyIntact.setupLogging(log_level: int, catch_sig: bool = True) -> None** Initialize logging Parameter ``log_level``: The max log level: 0 is info, -1 is warning, -2 is errors Parameter ``catch_sig``: Whether Intact should gracefully handle signals. For testing purposes where exit code of the program is important, set to false.