Homat
\renewcommand{\vec}[1]{\boldsymbol{#1}}
Introduction¶
Multiscale analysis of components with a heterogeneous microstructure¶
Material properties characterize the functionality of almost any product. For many industrial applications the final material properties result from a number of process steps. Each of them allows the tailoring of the microstructure in order to achieve the desired functionality. Bridging the scales is an essential feature of the multiscale analysis of components with locally different material properties ^{1}. In the case of multiscale material modeling and in homogenization in general, one proceeds from the detailed description of all material heterogeneities at the lower scale upward in order to obtain effective material properties accurately. Alternatively, in process modeling it is important to be able to step down to the scale of the real, heterogeneous material microstructure. The corresponding technique is known as the localization method. This coupling between the scales corresponds to a serial, handshake type exchange of data (e.g. stresses, strains, etc.) passing up and down the scales.
In Fig. 1 the multiscale analysis approach is outlined, using the example of the UO forming production step of an X65 steel pipeline tube with a pearliteferrite microstructure. The corresponding effective thermal properties (thermal conductivity, capacity, density) and thermoelastic properties (Hooke matrix, Young and shear moduli, Poisson ratios, thermal expansion coefficient, etc.) are evaluated by the homogenization tool HOMAT by applying the asymptotic homogenization method to heterogeneous materials. The theoretical background of this method is briefly summarized in Appendix Homogenization of Materials with Periodic Microstructure.
The effective nonlinear materials properties, like the elastoplastic flow curves, are determined via uniaxial virtual testing by applying, at the moment, the commercial finite element software package ABAQUS.
Then, the simulation of the UO forming process of the pipeline tube is performed at the macroscale and critical areas are detected. Locally in such a critical area, the macrostrain states are extracted and applied as boundary condition to the ferritepearlite microstructure in order to determine the corresponding distribution of microstrains and stresses on the Representative Volume Element (RVE) and the location of possible failure. The governing equations of the localization procedure within in the framework of the asymptotic homogenization method, are briefly outlined in Appendix Postprocessing of macroscale results: the localization step.
The homogenization and virtual testing tools are based on an accurate geometrical description of the complex microstructure (i.e. the topology of the individual phases and their spatial distribution) and the properties of the individual phases. The strength of these methods is their ability to consider transverse isotropic, cubic, orthotropic or generally anisotropic material properties. The phase properties may be either given by abinitio calculations or provided by models, experiments and/or databases. The geometry of the microstructure can be either simulated directly by the multicomponent phasefield program MICRESS or by the 3D cellular automata SphaeroSim for semicrystalline thermoplastics ^{2}^{,}^{3}^{,}^{4}, or specified experimentally via 3D CT scans or SEM micrographs ^{5}. As mentioned in Appendix Homogenization of Materials with Periodic Microstructure, not the entire microstructure is discretized at the microscale but only a RVE. For a material with periodic microstructure its size is defined by the period while for materials with random microstructure 1D spectral analysis or 2point probability functions are used ^{5} to specify a first guess of the RVE size. A more detailed description about the adopted stochastic homogenization procedure for materials with random microstructure is given in the bookchapter "Prediction of effective properties" ^{1}.
The homogenization method further requires the automatic generation of interface elements with zero thickness at the boundaries between phases, grains, spherulites or material domains. These interface elements allow the calculation of implicit forces/fluxes induced by the discontinuity of the material properties at these boundaries. To achieve this task a dedicated interface program, named Mesh2HOMAT, is supplied (see Fig. 2).
The preprocessor tools Mesh2HOMAT and Mesh2ABAQUS¶
The interface program Mesh2XXX has a generic name and it includes both specific interface tools: Mesh2HOMAT and Mesh2ABAQUS. The aim of this interface program is twofold. Firstly, it allows the conversion of a structured or unstructured mesh from one data format to another. Four different input file formats are currently available: the standard platform exchange format VTK, the universal format UNV of IDEAS, the Abaqus input format INP ^{6} and the MICRESS format, for which the possibility to read the compressed grain information directly is retained.
This preprocessor automatically detects the interfaces between adjacent material domains, grains or spherulites and generates new nodes and the corresponding interface elements. It evaluates the outward normal to the free surfaces and defines the pair of periodic surfaces as a function of the user specifications. Depending on the actual problem, either 1D, 2D or, by default, 3D periodic boundary conditions can be specified on the unit cell. The microstructure solutions of MICRESS and SphaeroSim, available in the initial VTK format, are rewritten in a VTK file for the newly generated, unstructured finite element mesh of the RVE. Eventually, Mesh2Homat writes a geometry file PROBLEM_NAME.hgeo
and a template of the HOMAT command file PROBLEM_NAME.hin
.
By specifying Mesh2ABAQUS as program name, this interface program generates automatically the six (or three in 2D) input files necessary to perform uniaxial virtual tests of the discretized RVE of the considered microstructure within ABAQUS. The program automatically generates the required periodic boundary conditions. It also reads the orientation of each grain of the simulated microstructure by MICRESS (PROBLEM_NAME.TabO
file) and transforms these orientations into a format readable by ABAQUS. The detailed description of the common and specific commands of both interface program versions Mesh2HOMAT and Mesh2ABAQUS is given in a separate manual.
Overview and Features of the Homogenization Tool HOMAT¶
Survey of HOMAT¶
The homogenization tool HOMAT is based on the asymptotic homogenization method applied to materials with periodic microstructures as outlined in Appendix Homogenization of Materials with Periodic Microstructure. This method allows the determination of the following effective properties:
 thermal conductivity, heat capacity and density
 Darcy permeability
 elastic stiffness (Hooke) matrix, thermal expansion coefficients, Young and shear moduli, Poisson ratios and volumetric eigenstrains
HOMAT started with the homogenization of cooled multilayer systems by the prediction of effective thermal properties of these systems ^{7}. Only periodic boundary conditions are specified on the RVE, even for materials with random microstructure. For these materials the periodicity is applied in a weak manner: nodes without a projection point in the same phase on the periodic face are eliminated from the periodic boundary condition. This is mostly relevant for Stokes flow homogenization ^{8}. The tool evaluates microscopic periodic displacements on the unit cell, which are used to express the unknown effective properties. Depending on the problem nature, up to 8 (thermoelasticity with eigenstrains) or 3 (thermal and fluid flow) boundary value problems are discretized by finite elements and solved on the unit cell successively. For the mechanical homogenization, the rigid body modes are fixed automatically and the user can specify the kind of the desired homogenization analysis: 3D, plane strain or plane stress. The righthand sides of these problems correspond to the material discontinuities at the grain or material boundaries plus the variation of the property inside each grain/domain (e.g. with the temperature or/and orientation of the domain). The evaluated effective properties of the unit cell are written in an ASCII file PROBLEM_NAME.hres
and in a material data file (.mat
) of ABAQUS type. Hence, the effective material properties can be directly transferred to the macroscale simulation programs (ABAQUS, DEFORM, SYSWELD, etc.). Moreover, for each righthand side and component, the distributions of the implicit fluxes, forces or boost forces and of the corresponding microscopic displacements on the unit cell are written in files of VTK format. In presence of mixed finite elements (modelling incompressible materials) also the periodic pressure fields are provided in VTK format. Since version 6.0, it is possible in thermoelastic homogenizations to write also the microscopic strain distribution on the RVE. In the case of Stokes flow homogenization, the periodic velocity and pressure fields are written by default. In order to reduce the storage data, these files are, by default, binary ones if the microstructure is simulated by MICRESS or SphaeroSim. These microscopic results can then be visualized with Paraview.
In the present user guide all the keywords of the input file PROBLEM_NAME.hin
and of the associated geometry file PROBLEM_NAME.hgeo
are detailed.
The program HOMAT is launched easily by typing:
homat PROBLEM_NAME
without any postfix on the command line or by dragging the command file PROBLEM_NAME.hin
on the icon/executable of HOMAT. Note that if you execute HOMAT without any model information on the command line, the program will display usage information:
Usage: homat_name <problem_name> <options>
<problem_name> : name of the homat input file without postfix
<options>
check : check of the consistency of the hin input file
The new option check
on the command line permits to check the consistency of the introduced input file PROBLEM_NAME.hin
. If activated, HOMAT will check if the required geometry and declared material properties and the input files (Microsim and/or Macro) are present. Moreover, it prints the chosen homogenization or localization analysis, the type of applied periodicity and some mesh informations, like the total number of nodes, the number of pure geometrical nodes and the number of additional interface nodes.
A detailed description how to use HOMAT in an efficient way is outlined and different homogenization and localization examples are provided with this version. Each example has its own definition file, allways named Def_Example.txt
. All the information required to run the considered example are outlined there. The specification of the command file PROBLEM_NAME.hin
for a discretized RVE of the considered heterogeneous material can be divided in different parts, which all have to be specified for a homogenization run:
 The definition of the thermal and/or thermoelastic material properties of each phase or domain of the microstructure have to be specified. A specific data file, named
MATERIAL.hmat
, has to be created for each phase or material. The material properties specific keywords will be described in Sec. Definition of the Materials Properties.  The geometry of the RVE is defined in the
PROBLEM_NAME.hgeo
file, created by Mesh2HOMAT. This file contains a list of the nodes, the bulk elements, the interface elements defined between phases or domains and the element faces belonging to outer surfaces. The geometry and model related keywords are presented in Sec. Definition of the Microstructure Model.  In Sec. Homat input files, the specific HOMAT keywords of the command file
PROBLEM_NAME.hin
are specified. This chapter involves also the description of the definition of localization runs in selected Gauss points of the preceding macrosimulation.  Finally, in Sec. HOMAT Output, in function of the specified homogenization analysis, the different kinds of output files are detailed and examples provided. The different solution files of a localization analysis are also explained.
Features of the HOMAT Version 6.0¶
In this section the key features of the HOMAT Version 6.0 are briefly outlined. The corresponding commands and input data are detailed in Sec. Definition of HOMAT specific Commands, dedicated to the specific HOMAT commands.
 HOMAT automatically evaluates the periodic displacement fields on the RVE (for example, 6 for a purely mechanical homogenization in 3D) required by the calculation of the desired effective properties.
 Several homogenizations can be defined and performed in one HOMAT run. For example, for a given microstructure RVE effective properties can be evaluated at several specified temperatures in the same HOMAT simulation.
 In the presence of a heterogeneous multilayer material, like the cooled combustion chamber wall composed of a superalloy, a bondcoat and a thermal barrier coating (see Fig. 25), the user can choose between multilayer or monolayer homogenization. In the case of monolayer homogenization, the user has to specify the layer for which effective properties are determined.
 Different kind of homogenizations (thermal, thermoelastic, ...) can be performed in the same HOMAT simulation. Note that the new Stokes flow homogenization cannot be combined with the two other kinds of homogenization.
 The periodic microscopic displacement fields \chi^{\text{rs}}, solutions of the linear problems [see Eq. \eqref{eq:solu_2}] and the implicit fluxes (thermal homogenization) and/or implicit forces (thermoelastic homogenization) are saved in files written in VTK format. In presence of mixed finite elements, used to accurately handle material incompressibility, the predicted microscopic pressure fields are also written in VTK format. In Stokes flow homogenization, the microscopic velocity and pressure fields are provided. The different solution files of HOMAT can be visualized by Paraview on the initial, undeformed unit cell or RVE.
 The considered periodic fields on the RVE (temperature, displacement, pressure, etc) are discretized in 3D by linear or quadratic finite elements of hexahedral, pentahedral and/or tetrahedral type: Since version 6.0, the initial element geometry can be also curved (using quadratic ansatz functions).
 In order to model accurately the quasiincompressibility of phases like the melt or the amorphous phase of a semicrystalline polymer, several mixed element types have been implemented in HOMAT 6.0: namely the equal order linear displacement and pressure elements (P1P1 and Q1Q1) stabilized by the double Gauss point integration, proposed by Li and He ^{9}^{,}^{2}, the mini element version, where the displacement field is augmented by a volume bubble. As a tetrahedron, this element discretization is automatically stabilized but not its hexahedral variant. Eventually, new stabilized elements with volume bubble and LiHe stabilization have been developed and tested successfully.
 The number of Degree of Freedoms (DOF) per node is now variable in the RVE model. If a quasiincompressible phase is present in combination with compressible ones, a 3D model will have nodes with 4 DOFs as well as other nodes with classical 3 DOFs.
 The elements can be either fully or reduced integrated ^{10}. By default, the reduced integration scheme is adopted. The user can specify the number of Gauss points per element to define the integration type.
 The predicted effective properties are written in a
PROBLEM_NAME.hres
file. This file serves for further postprocessing the homogenization results and also to print them in a convenient way. Since version 6.0, the user can write the effective properties also in a specified local axis system, e.g. In the grain axis system. Moreover, the effective material properties are written in aPROBLEM_NAME.mat
file of ABAQUS type. Thus, these predicted material properties can be directly included in the finite element model of the specimen at the macroscale.  The mechanical phase or domain properties can be either isotropic, transverse isotropic, cubic, orthotropic, monoclinic or anisotropic. In thermoelasticity, the thermal expansion can be either isotropic or orthotropic; whereas the eigenstrain is assumed to be isotropic. The thermal phase and domain properties can be now isotropic, transverse isotropic, orthotropic or anisotropic.
 If the phases or domains of the considered RVE are cubic with cubic axes corresponding to the RVE axis system, the effective Young and shear moduli, the effective Poisson and Zener anisotropy coefficients of the heterogeneous materials are expressed and written in the
PROBLEM_NAME.hres
result file.  Dependent on the nature of the homogenization problem, the user can choose in thermoelasticity to perform either a 3D homogenization (default case) or a 2D one of plane strain or stress type. In the case of thermal or fluid flow homogenization, either 3D or 2D homogenization can be selected.
 On the RVE 3D periodic boundary conditions (B.C.) are, by default, specified. Dependent on the nature of the heterogeneous materials of the multilayer system, 2D or 1D periodic B.C. can be also specified.
 HOMAT automatically detects in the list of the periodic nodes, the double nodes on wedges, which built the intersection of periodic surfaces, and removes the rigid body modes. The user has no specific fixations to define! Starting from a node on one wedge, a special link procedure is developed for 2D and 3D periodicity, like shown in Fig. 3. As the finite element mesh on one surface of a periodic pair can be different from the mesh on the other surface, we obtain a chain of linked D.O.F, which must be open in order to ensure equal movement of the linked nodes. In presence of 2D periodicity the corner nodes of the lowest face are fixed, like shown in Fig. 3 (left), whereas in the case of 3D periodicity all 8 corner nodes are fixed and three chains of linked wedge nodes are build (Fig. 3, right).
 At the interface between two material domains, tied contact is achieved. The slave nodes are condensed, only the master nodes are retained in the linear system. Thus, the system to be solved is reduced. This procedures leads to save computational time. After resolution, the results are duplicated on the condensed slave nodes.
 In the version 6.0 the concept of masterslave nodes is extended to the periodic surfaces Independent, if periodic surfaces have directly corresponding nodes or not, the slave nodes are eliminated and replaced, via a projection algorithm, detailed in Appendix Periodic boundary conditions via a masterslave projection algorithm, either by the directly corresponding master node or by the master nodes of the projection element. This way, the total system to be solved is reduced further.
 A major difference between the version 6.0 and the previous versions 5.x is the new treatment of the periodic B.C. via a masterslave projection algorithm, which doesn't use any penaty factor. This way, no additional stiffness at periodic surface nodes is introduced leading to larger microscopic displacement variations over the RVE and, subsequently, to larger corrections of volume averaging material property (see Eqs. \eqref{eq:eff_hoke_mat} and \eqref{eq:eff_eigenstress}). Therefore, the effective material properties predicted with HOMAT version 6.0 differ from previous predictions (version 5.3).
 The CSR sparse matrix format is now used to store the global system matrix. In this format, flexible number of DOFs per node are more easily treated and the implementation, in near future, of direct solvers is facilitated.
 An improved iterative BiCGstab solver is used to solve the linear systems of the mechanical, thermal or flow problem. The solver is enhanced by an automatic restart option and a new symmetric preconditioning, based on a scaling of the diagonal terms. This scaling is mainly important for the calculation of the displacement field \psi (calculation of the effective thermal expansion coefficients).
 The application of the periodicity B.C. on a heterogeneous material with random microstructure, like an opencell metallic foam, is delicate. For each periodic pair, the nodes belonging on the same phase or material domain on both surfaces are linked classically. But, nodes not having a projection point/element in the same phase on the master periodic surface are eliminated from the periodic boundary condition. This way, weak periodic B.C. are specified. When evaluating of the Darcy permeability, a special box can be defined inside the RVE where the velocities are used in the averaging process ^{11}. Thus, the main affected boundary layers near the free periodic surfaces can be omitted in the effective permeability evaluation.
Homat input files¶
Common features of HOMAT input files¶
A HOMAT input file is a plain ASCII data file. A template for the parameters of the keywords of the command file PROBLEM_NAME.hin
is created by the preprocessor Mesh2HOMAT. The user has to verify if the default options are adequate for the considered homogenization or localization run. He can modify and adapt the predefined parameters by using any text editor. The input files consist of a series of lines containing keywords, parameters and data. The following rules are applicable to all keywords and data lines of the command, geometry and material files read by HOMAT:
 The first nonblank character of each keyword line must be the slash symbol "
/
".  The keyword must be followed by a comma "
,
" if any parameter is specified.  Parameters must be separated by commas.
 Blanks on any line are ignored.
 A line is always limited to 256 characters, including blanks.
 Any line beginning with the hash symbol "
#
" is a comment line. Such lines can be placed anywhere in the input files of HOMAT, they are ignored by the “read line interpreter”.  Keywords and parameters are not case sensitive and be either written in upper or lower case letters, except file names.
 If a parameter has a value, the equal sign "
=
" is used. The value can be either an integer, a float number or a character string, depending on the parameter's meaning.  Continuation of a keyword line is sometimes necessary: for example, in presence of a large number of parameters of the considered keyword. If the last character on a keyword line is a comma, the next line is interpreted automatically as a continuation line of the previous line.
 Character strings can be built of up to 80 characters and are not case sensitive.
 Labels such as group or surface names are case insensitive. Their unique name has max. 80 characters and must not contain blanks.
 The name of the keywords and their parameters have been chosen as similar as possible to the keywords of the ABAQUS input file
.inp
. It is clear that HOMAT has its own keywords but mainly the keywords describing the material and geometry features of the model are made compatible to the ABAQUS input file. An easy change of the slash "/
" into "*
" allows an easy transformation of a HOMAT model in an ABAQUS model.
Definition of the Materials Properties¶
In this chapter the thermal and thermoelastic properties per phase or per domain are presented respectively. Per property class they are listed in alphabetical order of the keywords. All these temperature dependent material properties are collected, per material, in a special file, named MATERIAL_NAME.hmat
. Note that properties of the mushy zone can depend on the fraction solid, f_{S}, as the temperature is in nonunique above T_{sol} due to recalescence phenomena ^{3}. This file can be included easily in the command file of a specific homogenization application via the /INCLUDE
command explained more in detail in Sec. Definition of HOMAT specific Commands. This inclusion happens after the definition of the material by the keyword /MATERIAL
and its parameter. This way, the same phase or constituent material properties are specified once for several homogenization runs.
/MATERIAL
, NAME =
name
, TYPE =
type
This keyword is part of the model definition of the heterogeneous material at the microscale level. It is used to indicate the start of a material definition. A material data block can be defined directly in the command file PROBLEM_NAME.hin
or imported via the /INCLUDE
command by loading an existing MATERIAL_NAME.hmat
file, involving all properties of the considered phase or material. A material data block is defined by all keywords between the /MATERIAL
line and either another /MATERIAL
line or a keyword line that does not specify material properties. If a property is specified more than once for the same material, always the last definition will be used. Example:
/MATERIAL, NAME=ferrite, TYPE=ISO
The /MATERIAL
keyword has two required parameters, the NAME
of the constituent material or phase, like ferrite, austenite, ... and the TYPE
of the specified material model. The name of the material must be defined and is an alias for the defined material properties. The parameter NAME
will be used in the /SOLID SECTION
options in order to assign the material to a group of elements (see the definition of these keywords in Sec. Definition of the Microstructure Model). The second parameter TYPE
is essential not only for thermoelasticity but, since version 5.3, also for the thermal homogenization. It specifies the material model used for the considered phase or material domain. Indeed, if the user specifies TYPE=ISO
, which means that the ferrite phase is considered in the homogenization run as an isotropic material, only the Young modulus and the Poisson coefficient are considered in this simulation even if the cubic elastic constants are also provided in the material properties file .hmat
or in the command file .hin
. To specify the bcc ferrite phase as a cubic symmetric phase, characterized by the elastic constants H_{11}, H_{12} and H_{44}, the user has to specify the material model ORTHO
. Considering the elastic cubic symmetric behavior of the ferrite phase will, of course, improve the homogenization results obtained with the isotropic model. Therefore, the correct choice of the appropriate model is essential and has a direct impact on accuracy of the homogenization results. It is important to note that, by specifying TYPE=ORTHO
for the ferrite phase, the present Young and shear modules and Poisson coefficients in the material file .hmat are ignored because the material is assumed to be of cubic symmetry in the considered homogenization run. Following material models can be specified via the TYPE
parameter:
TYPE=ISO
: isotropic material behavior, defined only by the Young modulus and the Poisson coefficient and a single value for the thermal conductivity;TYPE=TR_ISO
: transverse isotropic material behavior, requiring the introduction of 5 independent material coefficients for the elastic behavior and two coefficients for thermal conductivity;TYPE=ORTHO
: orthotropic material behavior, requiring usually the definition of 9 independent material quantities and three values (or two in 2D homogenization) for the thermal conductivity;TYPE=ANISO
: anisotropic material behavior, requiring the definition of 21 material coefficients and 6 values for the thermal conductivity tensor;TYPE=SPHERU
: spherulite model for semicrystalline thermoplastics. This model is detailed in the recent publication ^{3}. Please contact ACCESS, if you want to use this special model in combination with SphaeroSim results.
Remark: The cubic symmetry is a subclass of the orthotropic material behavior and requires only the introduction of three material constants: User friendly introduction of this material behavior can be achieved via the /ELASTIC
keyword (see Sec. Thermoelastic properties). The specification of material models with mechanical anisotropy but thermal isotropy is delicate. In a simulation including thermal and thermoelastic homogenization, the material has to be specified either TR_ISO
, ORTHO
or ANISO
. Subsequently the definition of the thermal conductivity has to be in agreement with the adopted material behavior and also be specified as TR_ISO
, ORTHO
or ANISO
phase. On the other side, if only a pure thermal homogenization will be performed, this phase or material can be specified as isotropic and only one thermal conductivity value is specified. In this case, only isotropic mechanical are accepted in the material file MATERIAL_NAME.hmat
. For example, the presence of an orthotropic or anisotropic Hooke matrix definition will lead to an error message. Therefore, it is recommended to specify all phase properties to the anisotropy level of its thermoelastic behavior.
Note: The thermal conductivity definition of some materials has been modified in some existing examples in order to respect the agreement with the specified material model. It concerns Nicrofer wires in the hollow_grid
example, which are assumed to be transverse isotropic, the Yttrium stabilized Zirconium thermal barrier of the S30F
example, which is orthotropic and the sapphire fiber, which is anisotropic.
Remark: The calculation of the effective Darcy permeability tensor uses only unitary viscosities and densities per fluid domain see ^{8}, the Stokes flow homogenization in HOMAT requires for each specified /SOLID SECTION
keyword the name of a defined material. Of course, the real mechanical or thermal anisotropy of the specified material is neglected in this analysis.
Thermal properties¶

/CONDUCTIVITY
,TYPE =
type
This keyword defines the heat conductivity of the considered phase or material. Since the HOMAT version 5.3, the mandatory parameter
type
specifies the nature of the thermal conductivity: isotropic, transverse isotropic, orthotropic or anisotropic. As already mentioned, the definition of the thermal conductivity has to be specified in agreement with the adopted material model (see keyword/MATERIAL
). All values of the thermal conductivity may be temperature or fraction solid dependent. As HOMAT is units independent, the user can specify any type of coherent unit system (SI, CGS). As millimeter is defined as default length, classically the unit for thermal conductivity is W/(mm K). Depending on the phase nature, following different thermal conductivity definitions are available via theTYPE
parameter:
TYPE=ISO
: isotropic thermal behavior
After the 1^{st} line with the keyword/CONDUCTIVITY
and the parameterTYPE
, follow at the next lines the isotropic thermal conductivity with the corresponding temperature or fraction solid. The temperatures or fraction solids must be introduced in an ascending manner. Below the 1^{st} temperature, here 0° C, the thermal conductivity is set to the 1^{st} value, here 0.009 W/(mm.K). Above T=1000°C (last value), it is set to 0.2205 W/(mm.K). Between two introduced temperatures, HOMAT interpolates linearly.Note: The above detailed definition of a temperature dependent material property is general and will not be repeated by each property except in special cases like the elastic constants.
# /CONDUCTIVITY, TYPE=ISO 0.009000, 0.0 0.009120, 20.0 0.010180, 100.0 0.011500, 200.0 ............ 0.022050, 1000.0 # # lambda in W/(mm * K) # /DENSITY 0.0087, 0.0 # # temperature in C # Rho in g/mm**3 # /HEAT CAPACITY, TYPE=ISO 0.003370, 20.0 0.003500, 200.0 0.003880, 400.0 ................ 0.007130, 1000.0 # # rhocp = density * specific heat in J/(mm**3 *K) # Data provided by John Fernihough (1997) #
TYPE=TR_ISO
: transverse isotropic thermal behavior
The transverse isotropic case is characterized by a plane of isotropy. Only in the direction perpendicular to this plane, the material presents a different behavior than in the plane. Classical examples are unidirectional fibers. In order to specify a transverse isotropic thermal behavior, the user has only to introduce two independent thermal conductivities: \lambda_{\text{t}}, \lambda_{\text{p}}, and the temperature T or fraction solid f_{\text{S}}.\begin{equation*} \lambda_{\text{t}}\,, \lambda_{\text{p}}\,, T \end{equation*}Note that the subscript t stands for “transverse” and subscript p for “plane”. Example 2 shows such a definition of transverse isotropic thermal conductivity for a Nicrofer wire.
/CONDUCTIVITY, TYPE=TR_ISO 0.0107, 0.0113, 20. 0.0117, 0.0127, 100. 0.0133, 0.0144, 200. 0.0149, 0.0160, 300. 0.0164, 0.0176, 400. 0.0180, 0.0192, 500. 0.0194, 0.0206, 600. 0.0209, 0.0222, 700. # # unit W/(mm*K) #
TYPE=ORTHO
: orthotropic thermal behavior
A material presents an orthotropic behavior if we can identify three orthogonal symmetry planes in a local axis system (e_{1}, e_{2}, e_{3}). In this local axis system, the thermal conductivity tensor of such a phase has only diagonal terms: \lambda_{11}, \lambda_{22}, \lambda_{33} and the offdiagonal terms are zero. To describe such a thermal material behavior, the user introduces:\begin{equation*} \lambda_{11}\,, \lambda_{22}\,, \lambda_{33}\,, T \end{equation*}The Example 3 describes the definition of the orthotropic thermal conductivity of an yttrium stabilized zirconium oxide layer ^{5}.
/CONDUCTIVITY, TYPE=ORTHO 0.001436, 0.001436, 0.001436, 20. 0.001397, 0.001397, 0.001397, 200. 0.001379, 0.001379, 0.001379, 600. 0.001359, 0.001359, 0.001359, 800. 0.001385, 0.001385, 0.001385, 1000. # # unit W/(mm*K) #
TYPE=ANISO
: anisotropic thermal behavior
In the full anisotropic case, the material has no specific symmetry plane, whereas in the monoclinic case, it has only one symmetry plane. The thermal conductivity tensor \lambda_{ij} of an anisotropic material in 3D is usually a full 3\times3 matrix. Due to its symmetry, this tensor has only six independent conductivity terms:\begin{equation} \label{eq:lambda_ten} \lambda = \ \begin{bmatrix} \lambda_{11} & \lambda_{12} & \lambda_{13} \\ \bullet & \lambda_{22} & \lambda_{23} \\ \bullet & \bullet & \lambda_{33} \\ \end{bmatrix} \end{equation}To introduce such an anisotropic thermal behavior of the material in a local axis system (e_{1}, e_{2}, e_{3}), the user has to specify the six conductivity terms plus eventually the temperature in following manner:
\begin{equation*} \lambda_{11}\,, \lambda_{22}\,, \lambda_{33}\,, \lambda_{12}\,, \lambda_{13}\,, \lambda_{23}\,, T \end{equation*}As example, the effective thermal behavior of pearlite, which is an eutectoid phase mixture and composed of alternate ferrite and cementite layers, is outlined for the 18CrNiMo76 steel. These effective anisotropic thermal conductivities are results of a thermal homogenization run at the nanoscale, where a bilamella of ferrite and cementite is homogenized. Please look to the example
ferr_pearl
and the reference ^{12} for more details./CONDUCTIVITY, TYPE=ANISO 3.29137E02, 5.54883E02, 4.96613E02, 4.8956E04, 2.4484E04, 3.88469E03, 27. ## # #unit W/(mm*K) #
Note that this definition of anisotropic thermal conductivity differs from the Abaqus definition.
Remarks: If the thermal behavior of the material is not isotropic, then the keyword
/ORIENTATION
has to be specified in order to define the directions of material orthotropy or anisotropy. Alternative keyword definitions that HOMAT interprets as/CONDUCTIVITY
are/THERMAL CONDUCTIVITY
and/HEAT CONDUCTIVITY
. 

/DENSITY
This keyword allows the definition of the mass density \rho of the considered material. Usually the mass density is specified in g/mm^{3}. This thermophysical property can be defined as temperature dependent:
\begin{equation*} \rho\,, T \end{equation*}In Ex. 1 of the superalloy CMSX4, a constant density (value = 0.0087 g/mm^{3}) is defined for all temperatures.

/HEAT CAPACITY
The user has two definition possibilities: either he specifies the heat capacity or the specific heat per unit mass of the material. Both physical quantities are linked by the simple relation:
\begin{equation} \text{heat capacity (rcp)} = \text{density (rho)} * \text{specific heat (cp)} \end{equation}In Ex. 1 of the superalloy CMSX4, the temperature dependence of the heat capacity is specified and defined in the reference unit system J/(mm^{3} K). In presence of phase changes it takes the corresponding latent heat into account and a strong variation of its value is observed. The user specifies
\begin{equation*} \text{rcp}\,, T \end{equation*}

/SPECIFIC HEAT
The specific heat allows an alternate definition of the heat capacity of the material. The specific heat per mass unit of the material is defined as:
\begin{equation} \label{eq:spec_heat} c_{p}\left( T \right) = \frac{\text{dU}(T)}{\text{dT}}\,, \end{equation}where U(T) is the internal energy density of the material per unit mass due to thermal effects only.
Remarks: The user must choose: either he specifies the heat capacity of the material or the specific heat but not both physical quantities for the same phase or constituent material at the same time. HOMAT evaluates the effective density by volume averaging on the RVE, heat capacity and specific heat of the heterogeneous material. The emissivity of the phase or material does not need to be specified. HOMAT neglects, the radiation exchange inside the discretized unit cell. The user specifies
\begin{equation*} c_{p}\,, T \end{equation*}
Thermoelastic properties¶

/ELASTIC
,TYPE =
type
This option is used to define the elastic properties of a phase or a constituent material. Its linear thermoelastic behaviour is defined by the generalized Hooke law, linking the stress tensor \sigma_{ij} to the total elastic strain tensor \varepsilon_{kl}:
\begin{equation} \label{eq:stress} \sigma_{ij} = H_{ijkl}(T)\left[ \varepsilon_{kl}(\vec{u})  \alpha_{kl}(T)\Delta T  \kappa_{kl} \right]\,. \end{equation}The Hooke matrix H_{ijkl}(T) is a fourth order elasticity tensor and the thermal expansion coefficient \alpha_{kl}(T) are both temperature or fraction solid dependent quantities. The volumetric eigenstrains \kappa_{kl} are isotropic, per phase, and generally present in polycrystalline materials. The keyword
/ELASTIC
outlines the definition of the Hooke matrix. Depending on the number of symmetry planes for the elastic properties, a phase or constituent material can be classified as isotropic (infinite number of symmetry planes) or anisotropic (no symmetry planes). Some materials have a restricted number of symmetry planes. For example, orthotropic materials have three orthogonal symmetry planes for the elastic properties and monoclinic materials only one symmetry plane. The number of independent components of the Hooke matrix depends on such symmetry properties. The level of anisotropy and the mode of defining the elastic properties are chosen by theTYPE
parameter, which must be always specified on the keyword line. If the material is not isotropic, then the/ORIENTATION
option has to be specified in order to define the directions of anisotropy.
TYPE=ISO
: isotropic elasticity
The elastic properties are completely defined by specifying only the Young's modulus, E, and the Poisson coefficient \nu. For each temperature^{26} T the user provides one line of data\begin{align*} E, \nu, T \end{align*}Indeed, the shear modulus G of an isotropic elastic material is given by
\begin{equation} \label{eq:shear_mod_iso} G = \frac{E}{2(1 + \nu)}\,. \end{equation}These two parameters can be specified as function of the temperature (see Example 5 of the ferrite phase).
Remark: If the elastic properties of a constituent material are specified isotropic then the thermal expansion coefficient must be also defined as isotropic.
/ELASTIC, TYPE=ISO 211670.0, 0.27, 20. 208220.0, 0.27, 95. .......... 137210.0, 0.27, 760. 95300.0, 0.27, 911. # # unit MPa # /EXPANSION, TYPE=ISO 1.28925E5, 20. 1.31538E5, 100. ................ 1.72930E5, 850. 1.74279E5, 911. #
TYPE=TR_ISO
: transverse isotropic case
A special class of orthotropic material behaviour is the transverse isotropic case, which is characterized by a plane of isotropy at every point of the material. Assuming the plane 12 is the plane of isotropy, transverse isotropy requires E_{11} = E_{22} = E_{\text{p}} and \nu_{12} =\nu_{\text{p}} and G_{\text{p}} = \frac{E_{\text{p}}}{2(1 + \nu_{\text{p}})}, where the subscript p stands for “inplane”; but also G_{13} = G_{23} = G_{\text{t}} and \nu_{31} = \nu_{32} = \nu_{\text{tp}} and \nu_{13} = \nu_{23} = \nu_{pt} and G_{13} = G_{23} = G_{\text{t}}, where subscript t stand for “transverse” to the plane of isotropy. Thus, while \nu_{\text{tp}} is the Poisson ratio that characterizes the strain in the plane of isotropy induces by a stress normal to it, \nu_{pt} characterizes the strain in the normal direction to the plane of isotropy resulting from a inplane stress. In general these quantities are not equal. They are related by:\begin{equation} \label{eq:nu_triso} \frac{\nu_{\text{tp}}}{E_{\text{t}}} = \ \frac{\nu_{\text{pt}}}{E_{\text{p}}} \end{equation}The corresponding compliance matrix \mathbb{S}, inverse of the Hooke matrix \mathbb{H}, has the following expression in the Voigt notation
\begin{equation} \mathbb{S} = \mathbb{H}^{1} = \begin{bmatrix} \frac{1}{E_{\text{p}}} &  \frac{\nu_{\text{p}}}{E_{\text{p}}} &  \frac{\nu_{\text{tp}}}{E_{\text{t}}} & 0 & 0 & 0 \\  \frac{\nu_{\text{p}}}{E_{\text{p}}} & \frac{1}{E_{\text{p}}} &  \frac{\nu_{\text{tp}}}{E_{\text{t}}} & 0 & 0 & 0 \\  \frac{\nu_{\text{tp}}}{E_{\text{t}}} &  \frac{\nu_{\text{tp}}}{E_{\text{t}}} & \frac{1}{E_{\text{t}}} & 0 & 0 & 0 \\ 0 & 0 & 0 & \frac{1}{G_{\text{p}}} & 0 & 0 \\ 0 & 0 & 0 & 0 & \frac{1}{G_{\text{t}}} & 0 \\ 0 & 0 & 0 & 0 & 0 & \frac{1}{G_{\text{t}}} \\ \end{bmatrix}\,. \end{equation}For each temperature^{26} the user introduces only five independent material properties
\begin{align*} E_{\text{t}}, E_{\text{p}}, \nu_{\text{tp}}, \nu_{\text{p}}, G_{\text{t}}, T \end{align*}to specify the transverse isotropic material behavior, like shown in Example 6 for a unidirectional reinforced glassepoxy composite.
/ELASTIC, TYPE=TR_ISO 43565.0, 12800., 0.29, 0.307, 4143.85, 20. # # unit MPa # /EXPANSION, TYPE=ORTHO 4.51477E6, 1.02361D05, 1.02361D05, 20. #
TYPE=ENG
: orthotropic elasticity specified by engineering constants
Linear elasticity of an orthotropic material is most easily defined by giving its engineering constants: three Young moduli E_{1}, E_{2}, E_{3}; the Poisson coefficients \nu_{12}, \nu_{13}, \nu_{23} and the three shear moduli G_{12}, G_{13} and G_{23} associated to the material's or phase's principal material directions. This set of 9 independent moduli defines following elastic compliance matrix\begin{equation} \mathbb{S} = \ \mathbb{H}^{\mathbf{ 1}} = \begin{bmatrix} \frac{1}{E_{1}} &  \frac{\nu_{21}}{E_{2}} &  \frac{\nu_{31}}{E_{3}} & 0 & 0 & 0 \\  \frac{\nu_{12}}{E_{1}} & \frac{1}{E_{2}} &  \frac{\nu_{32}}{E_{3}} & 0 & 0 & 0 \\  \frac{\nu_{13}}{E_{1}} &  \frac{\nu_{23}}{E_{2}} & \frac{1}{E_{3}} & 0 & 0 & 0 \\ 0 & 0 & 0 & \frac{1}{G_{12}} & 0 & 0 \\ 0 & 0 & 0 & 0 & \frac{1}{G_{13}} & 0 \\ 0 & 0 & 0 & 0 & 0 & \frac{1}{G_{23}} \\ \end{bmatrix}\,. \end{equation}The Poisson coefficient \nu_{ij} measures the transverse strain in j direction, when the material is stressed in i direction. The symmetry of the Hooke and compliance matrices implies
\begin{equation} \label{eq:poisson_rat} \frac{\nu_{\text{ij}}}{E_{i}} = \frac{\nu_{\text{ji}}}{E_{j}}\,. \end{equation}Of course, the engineering constants can also be given as functions of the temperature, if necessary. For each temperature^{26}, the user provides two lines
\begin{align*} &E_{1}, E_{2}, E_{3}, \nu_{12}, \nu_{13}, \nu_{23}, G_{12}, G_{13},\\ &G_{23}, T \end{align*}as shown in Example 7 for an yttrium stabilized zirconium oxide layer ^{5}.
# /ELASTIC, TYPE=ENG 5.096191E+04, 3.124513E+04, 6.37023875E+04, 2.705990E01, 3.442350E01, 3.434799E01, 1.386069E+04, 1.386069E+04, 1.386069E+04, 700. 5.077967E+04, 3.121674E+04, 6.34745875E+04, 2.698298E01, 3.424140E01, 3.416010E01, 1.374312E+04, 1.374312E+04, 1.374312E+04, 750. .............. 5.011533E+04, 3.099163E+04, 6.26441625E+04, 2.682793E01, 3.369597E01, 3.351691E01, 1.351601E+04, 1.351601E+04, 1.351601E+04, 900. 4.975417E+04, 3.087421E+04, 6.21927125e+04, 2.673659E01, 3.347734E01, 3.331901E01, 1.349206E+04, 1.349206E+04, 1.349206E+04, 950. # /EXPANSION, TYPE=ORTHO 9.579838E06, 8.633990E06, 8.215011E08, 700 9.278161E06, 9.050553E06, 5.160571E07, 750 .................................. 9.449730E06, 9.183234E06, 3.867535E07, 900 9.380581E06, 9.157542E06, 5.460644E07, 950 #
As the elastic behaviour is orthotropic, the expansion coefficient of the material must be also specified as orthotropic.
Note: An alternative and totally equivalent definition of the elastic orthotropic material behaviour can be achieved via the combined definition of orthotropic Young and shear moduli and of orthotropic Poisson coefficients (please see the corresponding keywords and Ex. 12).

TYPE=ORTHO
: orthotropic elasticity via terms of the elastic Hooke matrix
Linear elasticity of an orthotropic material can also be defined by giving directly the nine independent elastic stiffness terms, as functions of the temperature, if desired. In this case, the Hooke matrix, written in Voigt notation, has following expression\begin{equation} \mathbb{H} = \begin{bmatrix} H_{11} & H_{12} & H_{13} & 0 & 0 & 0 \\ H_{21} & H_{22} & H_{23} & 0 & 0 & 0 \\ H_{31} & H_{32} & H_{33} & 0 & 0 & 0 \\ 0 & 0 & 0 & H_{44} & 0 & 0 \\ 0 & 0 & 0 & 0 & H_{55} & 0 \\ 0 & 0 & 0 & 0 & 0 & H_{66} \\ \end{bmatrix}\,. \end{equation}The symmetry of the Hooke matrix \mathbb{H} implies H_{ij} = H_{ji}. In text books about composite and anisotropic materials ^{13}^{,}^{14}, the reader can find not only the expressions linking the stiffness terms of expression (8) to the engineering constants (E_{i}, G_{ij} and \nu_{ij}) but also the stability criteria which these stiffness terms must fulfill in order to represent a material with positive deformation energy. The user has to specify two lines for each temperature^{26}:
\begin{align*} &H_{11}, H_{12}, H_{22}, H_{13}, H_{23}, H_{33}, H_{44}, H_{55},\\ &H_{66}, T \end{align*}The orthotropic elastic properties of the cementite phase of X65 pipeline steel are given as Ex. 8.
/ELASTIC, TYPE=ORTHO 388000., 156000., 345000., 164000., 162000., 322000., 15000., 134000., 134000., 23. # # Reference: Jiang et al. Jnl. of applied Physics, vol. 103, 043502, 2008. # /EXPANSION, TYPE=ORTHO 0.000006, 0.000006, 0.000006, 23. 0.000007, 0.000007, 0.000007, 200. ........ 0.0000190, 0.0000190, 0.000190, 1000. #

TYPE=CUBIC
: cubic symmetric elasticity
Metal phases, like the bcc ferrite or the fcc austenite in steels, present often with cubic symmetry, which is a special subclass of orthotropic materials. In this case, the diagonal stiffness terms H_{ii} (i=1,3) are identical; the shear stiffness terms H_{ii}, i=4,6 are also identical as well as the offdiagonal terms H_{ij}, i\neq j, i,j=1,3. In fact, the material is defined only by 3 independent elastic constants: H_{11}, H_{12} and H_{44}. The user has only to specify:\begin{equation*} H_{11}, H_{12}, H_{44}, T \end{equation*}/ELASTIC, TYPE=CUBIC 233100., 135440., 117830., 27. 228425., 131744., 115467., 100. .............. 127525., 73550., 64463., 900. # # Ref.: Ghosh & Olson, Acta Mat. 50 (2002) 265575 plus experimental data from # Rayne & Chandrasekhar, Phys. rev, 122 (1961) 171416. #
TYPE=ANISO
: fully anisotropic elasticity
For a fully anisotropic elastic material, 21 independent elastic stiffness parameters have to be specified. The symmetric Hooke matrix, in Voigt notation, of an anisotropic material is given by:\begin{equation} \label{eq:H_aniso} \mathbb{H} = \begin{bmatrix} H_{11} & H_{12} & H_{13} & H_{14} & H_{15} & H_{16} \\ & H_{22} & H_{23} & H_{24} & H_{25} & H_{26} \\ & & H_{33} & H_{34} & H_{35} & H_{36} \\ & S & & H_{44} & H_{45} & H_{46} \\ & & Y & & H_{55} & H_{56} \\ & & & M & & H_{66} \\ \end{bmatrix} \end{equation}In order to define this Hooke matrix, the user has to specify the 21 stiffness terms for each temperature^{26} in following manner:
\begin{align*} &H_{11}, H_{12}, H_{22}, H_{13}, H_{23}, H_{33}, H_{14}, H_{24},\\ &H_{34}, H_{44}, H_{15}, H_{25}, H_{35}, H_{45}, H_{55}, H_{16},\\ &H_{26}, H_{36}, H_{46}, H_{56}, H_{66}, T \end{align*}Note: HOMAT makes internally a further differentiation: it distinguishes between fully anisotropic material behavior with 21 independent Hooke matrix terms and the monoclinic case, with a symmetry plane either parallel to the Z=0 plane or to the X or Y planes. In these special cases, only 13 independent material coefficients have to be specified. In the thermoelastic case, the orthotropic thermal expansion coefficients \alpha_{i}\,,\text{ for } i=1, 3 must be specified. A full anisotropic thermal expansion tensor \alpha_{ij} cannot be introduced at the moment. As example, the monoclinic elastic behaviour of a monocrystalline sapphire fibre is outlined. The present isotropic Young modulus and Poisson coefficient are ignored when the parameter
TYPE
of the/MATERIAL
option is set toANISO
; whereas the/ELASTIC
keyword is ignored ifTYPE
is set toISO
./YOUNG MODULUS, TYPE=ISO 414000.0, 20.0 407000.0, 300.0  390000.0 900.0 385000.0, 1200.0 # # unit 1 N/mm2 = 1 Mpa # /POISSON RATIO, TYPE=ISO 0.27, 0.0 # /ELASTIC, TYPE=ANISO 456900., 149400., 456900., 106100., 106100., 462200., 0.0, 0.0, 0.0, 153700., 0.0, 0.0, 0.0, 20700., 135400., 20700.0, 20700., 0.0, 0.0, 0.0, 135400., 1100. # # thermal expansion coefficient 1/(Grad Celsius) function of temperature /EXPANSION, TYPE=ORTHO 0.0000069, 0.0000033, 0.0000033, 100.0 0.0000083, 0.0000077, 0.0000077, 500.0 0.0000090, 0.0000083, 0.0000083, 1000.0 0.0000095, 0.0000088, 0.0000088, 1500.0 #


/EXPANSION
,TYPE =
type
This keyword defines the thermal expansion coefficients of a phase or constituent material of the RVE. They are interpreted as total expansion coefficients with respect to a reference temperature T_{ref}, i.e. the thermal strain \epsilon_{th} of a material at final temperature^{26} T is given by
\begin{equation} \label{eq:th_expasion} \varepsilon_{\text{th}} = \mathbf{\alpha}(T)(T  T_{\text{ref}})\,, \end{equation}where \alpha(T) is the thermal expansion tensor at temperature T. As the objective of HOMAT is to determine the effective thermal expansion of the heterogeneous material, the temperature difference is always set to 1°C without loss of generality. Therefore, we have: \varepsilon_{\text{th}} = \alpha(T) The optional parameter
TYPE
on the keyword line specifies the adopted material model:ISO
,TR_ISO
orORTHO
. By default,TYPE=ISO
.TYPE=TR_ISO
has to be specifies for the transverse isotropic elastic behaviour, whereasTYPE=ORTHO
has to be defined for an orthotropic, cubic symmetric, fully anisotropic or monoclinic elastic material behaviour. Depending on theTYPE
value, the user specifies for each temperature^{26}:
TYPE=ISO
\begin{equation*} \alpha, T \end{equation*} 
TYPE =TR_ISO
\begin{equation*} \alpha_{\text{t}}, \alpha_{\text{p}}, T \end{equation*} 
TYPE=ORTHO
\begin{equation*} \alpha_{11}, \alpha_{22}, \alpha_{33}, T \end{equation*}
where \alpha_{11}, \alpha_{22} and \alpha_{33} are the thermal expansion in the directions 1, 2 or 3 of the local material axis system of the material domain (specified by the
/ORIENTATION
option) respectively; for the transverse isotropic case, \alpha_{\text{t}} is the thermal expansion coefficient in the transverse direction and \alpha_{\text{p}} the thermal expansion coefficient in the isotropic plane. The unit of the thermal expansion coeffcient is usually °C^{1}. Several definition examples have been already provided (see Ex. 58 and 10).Note: Alternatively to
/EXPANSION
the user can introduces the keyword/THERMAL EXPANSION
. Only the diagonal terms of the thermal expansion tensor can be introduced at the moment. Offdiagonal terms are neglected. 

/MOLAR VOLUME
This keyword defines the volume occupied by a mole of the considered phase/material. As mm is the default geometry dimension, its unit is mm^{3}/mol. In presence of different phases, this physical quantity must be specified for each phase in order to evaluate the corresponding isotropic, volumetric eigenstrain. For example, it occurs during a solid or liquid phase transformation, like the austenite to ferrite one in heat treatment processes. This quantity is directly linked to the phase/material density and to the molar weight and is, of course, temperature^{26} dependent. No additional parameter has to be specified. The user enters simply:
\begin{equation*} V_{mol}, T \end{equation*}The temperature dependent molar volume of the ferrite phase is given in Ex. 11.
# /MOLAR VOLUME 7099.52, 20. 7121.51, 100. ............ 7368.35, 850. 7391.90, 911. # expressed via the density and the molar weight of the alloy: 55.8309 # in mm^3/mol #

/POISSON RATIO
,TYPE =
type
This keyword is used in the case of isotropic, transverse isotropic or orthotropic material behaviour, defined by engineering constants, to introduce the Poisson coefficients of the considered phase or constituent material. Again in accord with the specified material model on the
/MATERIAL
keyword, the parameterTYPE
have to be defined on the keyword line. Classically, the user specifies:
TYPE=ISO
\begin{equation*} \nu, T \end{equation*} 
TYPE=TR_ISO
\begin{equation*} \nu_{\text{tp}}, \nu_{\text{p}}, T \end{equation*} 
TYPE=ORTHO
\begin{equation*} \nu_{12}, \nu_{13}, \nu_{23}, T \end{equation*}
It is important to note that, for a correct isotropic material behaviour definition, the option
/POISSON RATIO
must be accompanied by the keyword/YOUNG MODULUS
, whereas in the case of transverse isotropic or orthotropic material behaviour, this option must be specified in conjunction with the keywords/SHEAR MODULUS
and/YOUNG MODULUS
. Note that for the isotropic transverse material behaviour, the definition of the transverse and inplane Poisson coefficients, \nu_{\text{tp}} and \nu_{\text{p}} respectively, are identical as for the keyword/ELASTIC
. Ex. 10 illustrates the definition of an “isotropic” sapphire fibre. The orthotropic case is outlined in Ex. 12 for the yttrium stabilized zirconia oxide used as thermal barrier coating.# Young modulus (T) MPa # /YOUNG MODULUS, TYPE=ORTHO 5.096191E+04, 3.124513E+04, 6.37023875E+04, 700. 5.077967E+04, 3.121674E+04, 6.34745875E+04, 750. ......................... 5.011533E+04, 3.099163E+04, 6.26441625E+04, 900. 4.975417E+04, 3.087421E+04, 6.21927125e+04, 950. # /POISSON RATIO, TYPE=ORTHO 2.705990E01, 3.442350E01, 3.434799E01, 700 2.698298E01, 3.424140E01, 3.416010E01, 750 .......................... 2.682793E01, 3.369597E01, 3.351691E01, 900 2.673659E01, 3.347734E01, 3.331901E01, 950 # # shear modulus (T) MPa # /SHEAR MODULUS, TYPE=ORTHO 1.386069E+04, 1.386069E+04, 1.386069E+04, 700. 1.374312E+04, 1.374312E+04, 1.374312E+04, 750. ....................... 1.351601E+04, 1.351601E+04, 1.351601E+04, 900. 1.349206E+04, 1.349206E+04, 1.349206E+04, 950. #


/SHEAR MODULUS
,TYPE =
type
The user needs to specify this option only in the case of transverse isotropic or orthotropic material behavior, defined by a combination of Young and shear moduli and the Poisson coefficients. The
TYPE
parameter equals eitherTR_ISO
orORTHO
. The shear modules are usually specified in MPa and can depend on the temperature^{26}. The user introduces these modules as follows:
TYPE=TR_ISO
\begin{equation*} G_{\text{t}}, 0.0, T \end{equation*}Note that a second value must be specified for technical reasons. It is best set to 0.0.

TYPE=ORTHO
\begin{equation*} G_{12}, G_{13}, G_{23}, f_{S} \end{equation*}Note: The definition of the shear planes, noted here 12, 13 and 23, are achieved via the
/ORIENRATION
option. This option must be specified for each orthotropic or anisotropic material/phase.
Ex. 12 illustrates also the definition of this keyword.


/YOUNG MODULUS
,TYPE =
type
The Young modulus describes the basic engineering concept of elasticity. The keyword
/YOUNG MODULUS
can be used to describe either an isotropic, transverse isotropic or orthotropic material behaviour, defined by engineering constants. The additional parameterTYPE
defines the nature of the material model. Three values are permitted:TYPE=ISO
(default),TYPE=TR_ISO
andTYPE=ORTHO
. Only if its value corresponds to the value specified in the/MATERIAL
option, then it will be used to specify the Hooke matrix of the material, otherwise it will be ignored. Note that this option has to be used in combination with the keyword/POISSON RATIO
for an isotropic material and with/POISSON RATIO
and/SHEAR MODULUS
for a transverse isotropic or orthotropic material behaviour. The Young moduli are usually specified in MPa and can depend on the temperature^{26}. The user introduces either:
TYPE=ISO
\begin{equation*} E, T \end{equation*} 
TYPE=TR_ISO
\begin{equation*} E_{\text{t}}, E_{\text{p}}, T \end{equation*} 
TYPE=ORTHO
\begin{equation*} E_{1}, E_{2}, E_{3}, T \end{equation*}where 1,2,3 specify the directions of orthotropic symmetry; E_{\text{t}} defines the transverse Young modulus and E_{\text{p}} the inplane Young modulus (see also
/ELASTIC
,TYPE=TR_ISO
).
Example 8 illustrates the definition of this keyword.

Definition of the Microstructure Model¶
In this chapter all the keywords describing the geometrical model of the RVE on the microscale are outlined in an alphabetical order. The geometrical model has been pretreated by Mesh2HOMAT in order to generate the needed interface elements at phase/material boundaries and the periodic boundary conditions. All the geometrical information about the microstructure is summarized in the PROBLEM_NAME.hgeo
file, which is read by HOMAT. This file can be included easily in the command file of a specific homogenization application via the /INCLUDE
command explained more in detail in Sec. Definition of HOMAT specific Commands. This inclusion happens directly after the header, the problem, output and unit definitions, as shown in Ex. 13. After the description of the geometrical information provided by Mesh2HOMAT, the specific HOMAT keywords related to the microstructure model, like the /SOLID SECTION
, /INTERFACE
or /PERIODICITY
options, are outlined in detail.
/HEADING
Date: 25.04. Time: 10:24:20
Total number of Elements: 296
Total number of Nodes: 732
Total number of new Nodes: 76
Total number of Materials: 2
#
/PROBLEM, TYPE=THELAS, HOM_MOD=0, PER_MOD=2
#
/OUTPUT, FORMAT=VTK, TYPE=ASCII, IOP=2, VOLAVRG=YES, ENG=NO, VISLAV=YES, LOC_AXIS=NO
U
#
/UNITS, GEOMETRY=1.0, LENGTH UNIT=MM, TEMP OFFSET=273.15, ANGLE=GRAD
#
/INCLUDE, NAME=ehz3D_iso.hgeo.gz
#
/MATERIAL, NAME=VERRE, TYPE=ISO
/INCLUDE, NAME=verre.hmat
#
/MATERIAL, NAME=EPOXY
/INCLUDE, NAME=epoxy.hmat
/SOLID SECTION, NUMBER=1, ELSET=FIBER, NAME=fibre, MATERIAL=VERRE, TEMP=NO, MECH=CALC, FLOW=NO
#
/SOLID SECTION, NUMBER=2, ELSET=RESIN, NAME=matrice, MATERIAL=EPOXY,
TEMP=NO, MECH=HOM, FLOW=NO
#
/INTERFACE, NAME=FIBER:RESIN, NUMBER=1,
SURFACE1=INTERFACE_FIBER:RESIN,
SURFACE2=INTERFACE_RESIN:FIBER
#
/PERIODICITY, NAME=XX+, TYPE=TRANS, NUMBER=1, START=X, TARGET=X+
3, 0.0, 0.0
#
/PERIODICITY, NAME=X+X, TYPE=TRANS, NUMBER=2, START=X+, TARGET=X
3, 0.0, 0.0
#
/PERIODICITY, NAME=YY+, TYPE=TRANS, NUMBER=3, START=Y, TARGET=Y+
0.0, 3, 0.0
#
/PERIODICITY, NAME=Y+Y, TYPE=TRANS, NUMBER=4, START=Y+, TARGET=Y
0.0, 3, 0.0
Geometrical model information provided by Mesh2HOMAT¶

/ELEMENT
,TYPE =
type
,ELSET =
elset
This keyword defines the volume finite elements (FE) in a specific RVE domain of the heterogeneous material. Two additional parameters must be specified:
TYPE
andELSET
. The parameterTYPE
specifies the kind of finite element used for the considered domain. The element geometries available in HOMAT are summarized in Tab. 1.type
description C3D4
4node linear tetrahedral element C3D6
6node linear pentahedral element C3D8
8node linear hexahedral element C3D10
10node quadratic tetrahedral element C3D15
15node quadratic pentahedral element C3D20
20node quadratic hexahedral element The user familiar with ABAQUS ^{6} has directly recognized that the same element nomenclature is used in HOMAT as in ABAQUS.
Note: The use of quadratic elements is recommended to special applications where the accuracy of the prediction of effective properties is mandatory or in the localization analysis, where accurate local stresses near material interfaces are needed. Their use will severely increase of the number of D.O.F and, subsequently, of the CPU time and the stored data. In order to model accurately the mechanical and fluid flow behaviour of quasiincompressible materials and fluids, like melt in solidification processes, amorphous phase of semicrystalline polymers, mixed linear volume finite elements have been implemented in HOMAT 6.0. For each volume element type (tetra, penta and hexahedron), following different mixed variants have been introduced, detailed here for tetrahedral elements:
C3D4E
4node equal order tetrahedral element with linear displacement and pressure fields. This element version is not stabilized. Therefore, its use is not recommended and restricted to comparison purpose in scientific and benchmark applications.C3D4L
4node mixed stabilized equal order tetrahedral element with bilinear displacement and pressure fields. Following Li & He ^{9}, this element version is stabilized by a double Gauss point integration of the pressure term, avoiding pressure oscillations.C3D4M
4node tetrahedral element where the linear displacement field is enriched by a 3^{d} order volume bubble term. Its pressure field is always linear. As Arnold et al. ^{15} has proven, this tetrahedral element, named mini element, satisfies the infsup conditions and is automatically stable. However, the extension of the mini element to hexahedral C3D8M or pentahedralC3D6M
elements is not straightforward, as, for example for the hexahedron, three internal nodes instead one for the tetrahedron has to be added in order to satisfy the infsup condition.C3D4B
4node mixed stabilized tetrahedral element with bilinear displacement field enriched by a volume bubble term and a bilinear pressure field. This “bubble” element version is further stabilized by the LiHe method, by applying a double Gauss point integration of the pressure term.
For elements with pentahedral geometry, we can specify following mixed element variants:
C3D6E
,C3D6L
,C3D6M
andC3DB
. For the hexahedrons, we can define the same variants, namely:C3D8E
,C3D8L
,C3D8M
andC3D8B
. Note that for penta and hexahedral geometries, only the LiHe and the bubble variants are stable. In Mesh2HOMAT the new optiondiscre
allows the user to specify the desired finite element discretization for each domain/phase. By default, the classical bilinear displacement elementsC3D4
,C3D6
andC3D8
are used. The numerical benchmark tests (homogenization of a silicate sphere in an epoxy matrix, permeability prediction of flow around a cylinder, ..) have shown that the stabilized bubLi variant (C3D4B
,C3D6B
andC3D8B
) with enriched displacement field, is the most accurate element version. It is recommended to use this element version in mechanical applications with quasiincompressible materials and in Stokes flow homogenizations.The parameter
ELSET
specifies the name of the group of elements, which will be introduced under the same keyword/ELEMENT
. They define a specific element set, which can be easily addressed by its unique name. In Ex. 14 it defines, for example, the glass fibre. The keyword/ELEMENT
can be repeated several times to describe for example a complex domain composed of different types of elements. Independent of the element type, an element of the considered set is defined by its identification numberID_ele
followed by the list of the element nodes. For a linear hexahedral element, we have:ID_el, N1, N2, N3, N4, N5, N6, N7, N8
In Ex. 14 a selection of the geometry file generated by Mesh2HOMAT of the simple RVE of a reinforced UD glassepoxy composite (see Fig. 4) is presented. The
/ELEMENT
keyword has been used twice: once to define the fibre and once to specify the matrix around the fibre.# /NODE, NSET=ALLNODES 1, 1.190000e+00 ,+0.000000e+00 ,+2.500000e+01 2, +1.190000e+00 ,+0.000000e+00 ,+2.500000e+01 3, 1.173767e+00 ,+1.958856e01 ,+2.500000e+01 ................... 731, +1.125508e+00 ,3.864357e01 ,+2.500000e+01 732, +1.125523e+00 ,3.863900e01 ,2.500000e+01 /ELEMENT, TYPE=C3D8, ELSET=FIBER 1, 39, 21, 1, 3, 131, 124, 112, 111 2, 20, 2, 32, 73, 94, 93, 113, 165 3, 49, 39, 3, 4, 141, 131, 111, 110 ........................................... 223, 421, 424, 423, 420, 499, 502, 501, 498 /ELEMENT, TYPE=C3D8, ELSET=RESIN 76, 209, 254, 189, 185, 267, 335, 289, 266 77, 249, 218, 186, 195, 330, 299, 282, 283 ............................................ 295, 586, 582, 583, 587, 647, 643, 644, 648 296, 594, 592, 593, 595, 655, 653, 654, 656 # /ELSET, NAME=ALLELEMENTS FIBER RESIN # /SURFACE, NAME=INTERFACE_FIBRE:RESIN 1 , 6 149 , 6 2 , 5 150 , 5 3 , 6 ......................... 163 , 3 164 , 5 165 , 5 166 , 5 167 , 5 # /SURFACE, NAME=INTERFACE_RESIN:FIBRE 81 , 5 229 , 5 80 , 6 ...................... 269 , 4 270 , 6 272 , 6 273 , 6 274 , 6 /SURFACE, NAME=X 78, 5 79, 5 ................. 253, 5 254, 5 /SURFACE, NAME=X+ 76, 6 77, 6 ................ 234, 6 235, 6 /SURFACE, NAME=Y 76, 4 78, 4 ............... 99, 4 100, 4 /SURFACE, NAME=Y+ 224, 3 226, 3 ............... 247, 3 248, 3

/ELSET
,NAME =
name
, [GENERATE
]This option is used to assign elements to a specific element set or to combine existing element sets to a new one. The parameter
NAME
must be specified on the keyword line. In the Ex. 14, theNAME
isALLELEMENTS
. An element set with the same name must not already exist. The set name is case insensitive. The user has three ways to specify anELSET
:
List of element identification numbers, maximum 16 entries per line:
/ELSET, NAME=CRIT_AREA 45, 56, 89, 95, 85, 206, 218

In order to facilitate the explicit definition of the element numbers, belonging to a specific set, an optional parameter
GENERATE
can be specified on the keyword line. If present, element ranges can be expressed by their initial element number, their final number and an increment./ELSET, NAME=E2, GENERATE 20, 60, 5
In this case, the elements 20, 25, 30, 35, 40, 45, 50, 55 and 60 have been assigned to element set
E2
. 
As shown in Ex. 14, an
ELSET
may be specified as the union of previously defined element sets. Each subset appears on a separate line. In the special case of Ex. 14, the new element setALLELEMENTS
is the union of the element setsFIBER
andMATRIX
.


/NODE
, [SCAL
]This keyword specifies the nodes with their 3D coordinates of a part of the RVE or of the whole RVE of a heterogeneous material. The coordinates are given in the structural axis system of the RVE. At the moment, no local axis system can be defined for the introduction of nodal coordinates. On the keyword line the optional parameter
SCAL
can be introduced. This parameter allows the multiplication of the coordinates by the value of this parameter. Note that this scale factor overrides the scale value introduced via parameterGEOMETRY
of the more general keyword/UNITS
(see Sec. Options introduced before the RVE description). The users specifies for example:/NODE, SCAL=0.001 1, 0., 0., .0. 2, 10., 0., 0. 3, 0., 10., 0.
The coordinates of nodes 13 are multiply by the factor 10^{3} in order to pass from a micrometre to a millimetre definition of the RVE geometry.

/NSET
,NAME =
name
, [GENERATE
]Analogous to the definition of the element set
/ELSET
, keyword/NSET
is used to assign nodes to a specific node set. The parameterNAME
, specifying the name of the considered node set, is required. A node set with the same name must not already exist. The set name is also case insensitive. The user has three ways to specify a specificNSET
:
List of node numbers, maximum 16 entries per line
/NSET, NAME=AROUND_HOLE 45, 56, 89, 95, 85, 206, 218

In order to facilitate the generation of the list of node numbers an optional parameter
GENERATE
can be specified on the keyword line. If present, nodes ranges can be expressed by their initial node number, their final number and an increment./NSET, NAME=FREE_SURF, GENERATE 12, 24, 2
In this case, the nodes 12, 14, 16, 18, 20, 22 and 24 have been assigned to the node set
FREE_SURF
. 
List of previously defined nodes sets  one per line  to specify a new node set. This way it is easy to join already defined node sets.


/SURFACE
,NAME =
name
This keyword is used to assign element faces to a specific surface group. The parameter
NAME
, specifying the name of the considered surface, is required. Surfaces are used to specify the interfaces between different phases, materials, where jumps of the thermophysical properties occur. They are also used to define the periodic boundary conditions at the free surfaces of the RVE and contact heat transfer between phases, constituents inside the RVE for thermal homogenization. The numbering of the faces, belonging to the surface, depends on the element type (see Figs. 5 and 6 and Tabs. 2  4).face nodes (linear) nodes (quadratic) 1 1, 2, 3 1, 2, 3, 5, 6, 7 2 1, 4, 2 1, 4, 2, 8, 9, 5 3 2, 4, 3 2, 4, 3, 9, 10, 6 4 3, 4, 1 3, 4, 1, 10, 8, 7 face nodes (linear) nodes (quadratic) 1 1, 2, 3, 4 1, 2, 3, 4, 9, 10, 11, 12 2 5, 8, 7, 6 5, 8, 7, 6, 18, 15, 14, 12 3 1, 5, 6, 2 1, 5, 6, 2, 17, 12, 18, 9 4 2, 6, 7, 3 2, 6, 7, 3, 18, 14, 19, 10 5 3, 7, 8, 4 3, 7, 8, 4, 19, 15, 20, 11 6 4, 8, 5, 1 1, 5, 8, 4, 17, 16, 20, 12 face nodes (linear) nodes (quadratic) 1 1, 2, 3 1, 2, 3, 7, 8, 9 2 4, 6, 5 4, 6, 5, 12, 11, 10 3 1, 4, 5, 2 1, 4, 5, 2, 13, 10, 14, 7 4 2, 5, 6, 3 2, 5, 6, 3, 14, 11, 15, 8 5 1, 3, 6, 4 1, 3, 6, 4, 9, 15, 12, 13 Independent of the element type, the faces are numbered in such a way that the normal to this face is always oriented inside the volume of the element. This face definition is completely identical to the notation used in ABAQUS. The user specifies on the keyword line the parameter
NAME
and its value. Then, on the following line the element number and the number of face is specified. This line is repeated until all faces belonging to the considered surface are introduced./SURFACE, NAME=FreeSurf_1 1, 6 2, 1
This entry assigns the nodes of face 6 of element 1 and the nodes belonging to face 1 of element 2 to the surface with the name
FreeSurf_1
. For the RVE of the UD reinforced composite of Fig. 4, the surfacesX
,X+
,Y
andY+
are specified on the free surfaces of the RVE, which are normal to the X or Y direction respectively. These surface definitions are required by the/PERIODICITY
keyword (see Sec. Model information defined in the command file). At the interface between the fibre and the matrix two surfaces namedINTEFACE_FIBRE:RESIN
andINTERFACE_RESIN:FIBRE
respectively, are automatically defined by Mesh2HOMAT. The nomenclature of the surface names of an interface between two grains, phases, or material regions follows always following rules: 1^{st} surface:
INTERFACE_{1st material region}:{2nd material region}
 2^{nd} surface:
INTERFACE_{2nd material region}:{1st material region}
where the definition of the 1^{st} and 2^{nd} material regions depends of their order in the list of ELSET's.
 1^{st} surface:

/XTUPLE
, *=MAX ==*/=nbr=/, *=NUM ==*/=nodes=/This new keyword provides for each geometrical master node, where several (> 2) materials are connected, the list of nodes with the same coordinates. The mandatory second parameter
MAX
specifies the maximum number of nodes having the same coordinates; the second parameterNUM
specifies the number of XTUPLE specified automatically by Mesh2HOMAT.For the
hollow_grid
example, Mesh2HOMAT writes in thegrid_mono_th.hgeo
file following/XTUPLE
information:# /XTUPLE, MAX=3, NUM=648 5, 32154, 32156 7, 32174, 32176 17, 32242, 32245 24, 32292, 32295 ..... 2176, 34294, 34295 2177, 34211, 34218 2178, 34258, 34262
Model information defined in the command file¶

/HEAT TRANSFER
,TYPE =
IDEAL
FCT
, [HVAL =
hval
]This keyword is only read in thermal homogenization simulations. It is a subkeyword related to a specific interface and must be introduced directly behind the concerned
/INTERFACE
. As millimeter is the default unit of length, the heat transfer coefficient is usually specified in W/(mm^{2} K). Dependent on the nature of the two regions in contact, the heat transfer coefficient can be assumed to be ideal or has a temperature dependent value, obtained from experiments. The user must define on the keyword line the nature of the heat transfer coefficient: eitherIDEAL
orFCT
(= function of the temperature).
Ideal heat transfer
/HEAT TRANSFER, TYPE=IDEAL, HVAL=0.008
In the case of ideal heat transfer between two contact regions, the user can defined a specific value corresponding to this ideal heat transfer. By default,
HVAL=0.01
W/(mm^{2} K). In the present example, it is defined a little lower: 0.008 W/(mm^{2}K). If the parameterHVAL
is specified, its value overrides the parameterHTR_IDEAL
, which is defined in the/HOM PARA
keyword (see Sec. Definition of HOMAT specific Commands) for the whole RVE. 
Function of the temperature
/HEAT TRANSFER, TYPE=FCT 0.0005, 500.0 0.0008, 700.0 0.001. 900.0
Linear interpolation with the temperature is achieved between these values. If the temperatures of both regions in contact differ, HOMAT takes, per convention, the temperature of the 1^{st} region in the region list as the reference temperature at which the interpolation happens.


/INTERFACE
,NAME =
name
,NUMBER =
number
,SURFACE1 =
surface1
,SURFACE2 =
surface2
This keyword is automatically generated by Mesh2HOMAT and allows the definition of an interface between two grains, phases or materials at the microscale. Both surfaces belonging to an interface with zero thickness and having double coincident nodes (one per surface in contact) must have a coherent surface mesh: the number of nodes and triangular or quadrangular surface elements must be exactly the same on both interface sides. The keyword line can be followed by a
/HEAT TRANSFER
option in the case of thermal homogenization. For the unique interface between fibre and matrix of Ex. 13, Mesh2HOMAT generates the following interface definition:/INTERFACE, NAME=FIBER:RESIN, NUMBER=1, SURFACE1=INTERFACE_FIBER:RESIN, SURFACE2=INTERFACE_RESIN:FIBER
On the keyword line, the parameter
NAME
is required and specifies the name of the interface. The name is always a combination of the names of the two domains in contact, hereFIBER
andRESIN
. Further a second parameterNUMBER
is introduced there. It specifies the identification number of the interface. In this simple example, this number equals 1. On the second and third line, the names of both surface groups, previously defined in thePROBLEM_NAME.hgeo
file, are introduced respectively. ParameterSURFACE1
specifies always the surface of the first specified region in the geometrical model description andSURFACE2
the 2^{nd} one. In Ex. 13 the fiber domain is the 1^{st} region and the resin the 2^{nd} region. Mesh2HOMAT generates such a specific interface for each combination of different regions, having their own properties. In the case of grains of a polycrystal, belonging to the same phase, e.g. ferrite or austenite, the material properties are the same, of course, but not their orientation. In this case also a specific interface is generated between adjacent grains. For complex polycrystal microstructures (> 100 grains) several interfaces (> 1000) have to be specified. In this case manual definition is hardly an alternative to automatic generation via Mesh2HOMAT.

/ORIENTATION
,NAME =
name
,DEFINITION =
definition
, [ROTX
=rotx
], [OUT_ORIE
]This option is used to define a local coordinate system for a region, grain in order to specify its transverse isotropic, orthotropic or anisotropic material properties or in order to specify the local axis system in which the calculated effective properties are output. On the keyword line the required parameter
NAME
assigns a name to each orientation definition. This name is referenced in the/SOLID SECTION
keyword (see below). The second parameter,DEFINITION
, specifies the nature of the introduced orientation definition. At the moment, HOMAT recognizes following orientation definitions:
DEFINITION=QUATERNION
A quaternion has four components and is composed of a vector plus a scalar: \mathbf{q} \triangleq \lbrack\mathbf{v}^{T},\ q_{4}\rbrack, whose norm is one: \left v \right^{2} + \ q_{4}^{2} = \ \mathbf{q}^{T}\mathbf{q} = 1. A quaternion representing without the rotation by an angle α around a unit vector \vec{a} in 3D has the following expression\begin{equation} \label{eq:quaternion} [a_{1}\sin(\alpha/2), a_{2}\sin(\alpha/2), a_{3}\sin(\alpha/2), \cos(\alpha/2)]\,. \end{equation}The user introduces for example:
/ORIENTATION, NAME=grain_12, DEFINITION=QUATERNION 0.880735, 0.2870, 0.33212, 0.177857

DEFINITION=EULER
In 3D any rotation of a structural axis system to a local axis system and back, can be described by a sequence of three Euler rotations about coordinate axes, where two successive rotations may not be done around the same axis. Several ways to define the Euler angles exist in the literature. HOMAT uses the convention shown in Fig. 7, where the first rotation is about the Z structural axis (angle \phi, rotation matrix \vec{D}); followed by a rotation about the new X' axis (angle \theta \in [0, \pi], rotation matrix \vec{C}) and finally a rotation about Z' is achieved (angle \psi, rotation matrix \vec{B}). Any kind of rotation matrix \vec{A} from a structural to a local axis system can then be expressed in 3D by following successive matrix products:\begin{equation} \label{eq:euler_ang} \vec{A} = \vec{B} \vec{C} \vec{D} \end{equation}For an Euler angle orientation definition in degrees, the user introduces following sequence:
/ORIENTATION, NAME=grain_8, DEFINITION=EULER 90, 14.0, 0.
where \phi=90°, \theta=14° and \psi=0°.

DEFINITION=ANGLE_AXIS
An orientation in the 3D space can also be defined by a unit axis, defined by its components in the structural system and by a rotation angle \alpha. In this case the user provides simply at first the angle in degrees followed by the three components of the rotation axis in the 3D space./ORIENTATION, NAME=grain_8, DEFINITION=ANGLE_AXIS 90.0, 0.0, 1.0, 0.0

DEFINITION=ROTMA
The user can also directly specify the unitary rotation matrix \vec{A} (eq. 12), which is composed of the scalar products of the unit vectors of structural (\vec{E}_{j}) and the new local \vec{e}_{i} axis system:\begin{equation*} a_{ij} = \vec{E}_{i} \cdot \vec{e}_{j} \end{equation*}The explicit expression of the rotation matrix \vec{A} is given by
\begin{equation} \label{eq:rot_mat} \vec{A} = \left( \begin{matrix} a_{11} & a_{12} & a_{13} \\ a_{21} & a_{22} & a_{23} \\ a_{31} & a_{32} & a_{33} \\ \end{matrix} \right)\,. \end{equation}The definition of a local orientation axis system via the rotation matrix \vec{A} is achieved by introducing its components. This way, a direct link to the output of EBSD and metallographic measurements is achieved. Note that, HOMAT performs the inverse transformation from the local axis system, in which the material properties are defined, to the structural one, in which the homogenization analysis is performed.
/ORIENTATION, NAME=grain_18, DEFINITION=ROTMA a11, a12, a13, a21, a22, a23, a31, a32, a33
or , for a particular example,
/ORIENTATION, NAME=grain_18, DEFINITION=ROTMA 0.29237, 0.0, 0.9563, 0.0, 1.0, 0.0, 0.9563, 0.0, 0.29237

DEFINITION=COORDINATES
Analogous to the orientation definition in ABAQUS, the user can also specify a local, rectangular coordinate axis system by defining the origin of this system \vec{c} and two points \vec{a} and \vec{b}, as shown in Fig. 8 Point \vec{a} must be on the local axis X' and point \vec{b} must be in the plane X'Y', not necessary located on axis Y'. Intuitively, the user selects point \vec{b} on or near the local axis Y' (see Fig. 8). If definitionCOORDINATES
is specified, points \vec{a} and \vec{b} are defined directly by introducing their coordinates. If the user wants to define an origin of the local axis different from the structural one, he must also specify the coordinates of a third point \vec{c}. Otherwise, the coordinates of the two nodes \vec{a} and \vec{b} suffice. In the general case with different origin, the user introduces following orientation definition sequence:/ORIENTATION, NAME=some_name, DEFINITION=COORDINATES a1, a2, a3, b1, b2, b3, c1, c2, c3
or , for a particular example,
/ORIENTATION, NAME=wire_hori, DEFINITION=COORDINATES 2.0, 0.0, 0.0, 1.0, 3.0, 0.0, 1.0, 0.0, 0.0

DEFINITION=NODES
Similar toDEFINITION=COORDINATES
, but the points \vec{a}, \vec{b} and (\vec{c}) are located at two (three) nodes of the finite element mesh of the analysed RVE. If orientation definitionNODES
is selected, the user specifies the nodes corresponding to points \vec{a}, \vec{b} and (\vec{c}) by their global node numbers in the mesh. The input sequence might look like this:/ORIENTATION, NAME=wire_vert, DEFINITION=NODES 12, 34, 54

DEFINITION=OFFSET TO NODES
In this special definition of a local axis system, points \vec{a}, \vec{b} and (\vec{c}) are defined by specifying the local node numbers, relative to an element. By default, the origin corresponds to the first node number of each element belonging to the considered element set. For example, local node number 4 is given as the location of point \vec{a} and local number 5 to define point \vec{b}. Then, the local X' axis is defined to be parallel to the (1, 4) side of the hexahedral element and local axis Y' is in the plane formed by the nodes 1, 4 and 5 of this element. In this case, the user introduces simply:/ORIENTATION, NAME=fabric, DEFINITION=OFFSET TO NODES 4, 5
Of course, if the origin of the local axis systems differs from the first element node, the user has to specify also point \vec{c} by the corresponding local origin node. This orientation definition is very flexible and powerful. Indeed, in presence of an orthotropic curved layer, this orientation definition allows the user to specify easily a different local axis system in each element of the considered element set. The curvature of the layer implies such a local orientation definition per element. Note that, as the local axis system varies from element to element, also the Hooke matrix, expressed in the structural axis system, varies in the considered material domain, leading to the generation of implicit body forces, which are taken into account in the homogenization analysis (see Eq. \eqref{eq:solu_2} of Appendix Homogenization of Materials with Periodic Microstructure).

ROTX=90
for 2D MICRESS microstructures
The additional optional parameterROTX
must be specified on the keyword line in the case of 2D homogenization of microstructure predicted by MICRESS and no specific orientation file (.TabO
file) is transferred to HOMAT. Indeed, in 2D phase field simulation MICRESS defines the discretization domain in the plane XZ. Structural axis Y is perpendicular to this plane. Whereas HOMAT specifies always a 2D plane strain or stress model to homogenize in the classical XY plane. Therefore, an additional rotation of \alpha = 90° about the common X axis must be achieved to turn the MICRESS model with its grain orientation specified in the \vec{E}_{X}, \vec{E}_{Z}, \vec{E}_{Y} axis system (see Fig. 9), in the e_{X}, e_{Y}, e_{Z} axis system of HOMAT. In 2D, the user must specify this additional rotation via theROTX
parameter:/ORIENTATION, NAME=grain_1, DEFINITION=QUATERNION, ROTX=90 0.880735, 0.287, 0.3321, 0.177857

OUT_ORIE
This additional optional parameter permits to define the local axis system in which the calculated effective properties are printed in the result files (PROBLEM_NAME.hres
andPROBLEM_NAME.mat
) in addition to the structural axis system of the RVE. In fact, it permits to save the name of the actual orientation definition for this additional output./ORIENTATION, NAME=or_grain3, DEFINITION=EULER, OUT_ORIE 45, 14.0, 60.
The orientation
or_grain3
, defined by the Euler angles \phi=45°, \theta=14° and \psi=60°, specifies in this example the local axis system for the output of the effective orthotropic or anisotropic material properties.Note: The
OUT_ORIE
parameter can be specified only once for a specific orientation of the whole homogenization model. The locally variable orientation definitionOFFSET TO NODES
cannot be used to specify the output local axis system.


/PERIODICITY
,NAME =
name
,TYPE =
TRANS
ROTA
,NUMBER =
number
,START =
start
,TARGET =
target
This important boundary condition must at least be defined in one direction, for example structural direction X, on the microstructure RVE. Classically the RVE is periodic in all three directions. As further outlined in Appendix Homogenization of Materials with Periodic Microstructure, the periodicity B.C. introduces a constraint between displacements belonging to two opposite surfaces separated by a distance called period. Mesh2HOMAT automatically identifies, for each structural direction the appropriate free surfaces. A special nomenclature is introduced: the surface normal to the X, Y or Z directions, whose X, Y or Z coordinates are the smallest, are named
X
,Y
orZ
; whereas the surfaces with the largest X, Y or Z coordinates are named theX+
,Y+
orZ+
surfaces. Depending on the adopted periodicity assumption, e. g. 3D, 2D or 1D, Mesh2HOMAT prepares the definition of the linked surface pairs. In Ex. 13 of the UD reinforced composite, 2D periodicity is specified, as the unit cell can be repeated in X and Y direction by a period equal to its extension, but not in the fibre direction, as the fibre is assumed to be endless. In this case, Mesh2HOMAT generates following:/PERIODICITY, NAME=XX+, TYPE=TRANS, NUMBER=1, START=X, TARGET=X+ 3.0, 0.0, 0. /PERIODICITY, NAME=X+X, TYPE=TRANS, NUMBER=2, START=X+, TARGET=X 3.0, 0.0, 0. /PERIODICITY, NAME=YY+, TYPE=TRANS, NUMBER=3, START=Y, TARGET=Y+ 0.0, 3.0, 0. /PERIODICITY, NAME=Y+Y, TYPE=TRANS, NUMBER=4, START=Y+, TARGET=Y 0.0, 3.0, 0.
After the keyword
/PERIODICITY
;, the parameterNAME
specifies the name attributed to the considered link between two periodic surfaces. This name is a composition of the names of the two linked surfaces. The first surface name corresponds here to the surface, which is projected on the other surface, specified by the second name. Then, the parameterTYPE
defines how the surfaces are linked. HOMAT allows the introduction of two different transfer modes: translation or rotation.
If
TYPE=TRANS
is specified, a translation vector, defining the period, is introduced on the second line. Its components specify, per direction, the increments need to go from the slave face (1^{st} surface) to the master surface on the opposite RVE side (2^{nd} surface). The required parameterSTART
specifies the name of the starting (slave) surface and the parameterTARGET
defines the master surface, on which the slave nodes are projected. 
If
TYPE=ROTA
is specified, the user provides in the XZ plane a rotation angle \theta_{y} (in degrees) around the common structural axis Y. This special feature allows to define periodic B.C. on curved cylindrical RVE's, is detailed in reference ^{16}./PERIODICITY, NAME=XX+, TYPE=ROTA, NUMBER=1, START=X, TARGET=X+ 0, 15.7., 0.
Note: at the moment, only a rotation \theta_{y} about the common axis Y is implemented in HOMAT. For each linked pair of periodic surfaces, for example X and X+, the program Mesh2HOMAT specifies always two
/PERIODICITY
commands: the 1^{st} one to define the translation/rotation in one direction and the other one to specify the translation/rotation in the opposite direction. Indeed, it is important to specify each surface one time as slave face and the other time as master face. As 2D periodicity is defined on the glassepoxy RVE of Fig. 4, two pairs of linked surfaces, in X and Y directions respectively, are generated with four/PERIODICITY
commands in Ex. 13.


/SOLID SECTION
,NUMBER =
number
,ELSET =
elset
,NAME =
name
,MATERIAL =
material
, [ORIENTATION =
orientation
, ]TEMP =
NO
HOM
CALC
,MECH =
NO
HOM
CALC
,FLOW =
NO
HOM
CALC
This keyword is used to define physical regions with common properties, defined previously by the
/ELSET
option and corresponding, for example, to grains of a polycrystal or to regions of the RVE of a heterogeneous material. Some parameters are always required and theORIENTATION
parameter must only be specified if the considered section is transverse isotropic, orthotropic or anisotropic (see for example keyword/ELASTIC
). On the command line the user specifies successively the following parameters:ELSET
: this parameter is set equal to the name of the element set containing all elements for which the material behaviour is being defined via this option. Note that alternatively, the parameterGROUP
can be specified. It has exactly the same meaning.NUMBER
: defines the actual position in the list of all specified solid sections. Its count starts from 1 and it is incremented at each new solid section definition.NAME
specifies the name of the considered solid section.MATERIAL
set this parameter equal to the name of the previously specified material to be applied to this solid section.
ORIENTATION
for the transverse isotropic, orthotropic and anisotropic materials the user has to specify the name of the previously defined local axis system via keyword/ORIENTATION
. Of course, this parameter can be neglected for isotropic materials.Remark: In order to reduce drastically the definition effort in the case of MICRESS RVE's with a lot of grains, having each another grain orientation, a direct import of the orientations from the
TabO
file can be specified via the parameterORIE
of option/MICROSIM
(see Sec. Definition of HOMAT specific Commands). In this special case the parameterORIENTATION
can be omitted in the/SOLID SECTION
definition.
In addition to these parameters common with ABAQUS, the following special HOMAT parameters have to be specified. Their definition is essential as they define if the considered solid section participates in the actual thermal and/or thermoelastic or Stokes flow homogenization run.
TEMP
: IfTEMP=CALC
orHOM
, then the considered solid section is part of the planed thermal homogenization analysis. IfTEMP=NO
, this solid section does not participate to the thermal homogenization analysis. This way, the simulation model can be easily reduced to solid sections of interest and omit non interesting regions of the RVE.MECH
: this parameter can have also the values:CALC
,HOM
orNO
. In the caseCALC
orHOM
, the considered solid section is part of the planed thermoelastic homogenization and, of course, in the caseNO
, it does not participate to the analysis.FLOW
: this parameter can have also the values:CALC
,HOM
orNO
. In the case ofCALC
orHOM
, the considered "solid section" is a fluid domain and participates to the Stokes flow homogenization. Please specify for all the solid domains,FLOW=NO
, as in these regions no flow happens.
Note: The new parameter
problem
of Mesh2HOMAT specifies not only the desired homogenization analysis but also provides, for each homogenization analysis type, a default definition of the parametersTEMP
,MECH
andFLOW
. Of course, the user can modify this definition to suit his needs. In Ex. 13 of the UD reinforced composite a thermoelastic homogenization analysis is defined. Therefore for both constituents, the fibre and matrix regions, a/SOLID SECTION
keyword is specified withTEMP=NO
andMECH=CALC
orHOM
:/SOLID SECTION, NUMBER=1, ELSET=FIBER, NAME=fibre, MATERIAL=verre, TEMP=NO, MECH=CALC # /SOLID SECTION, NUMBER=2, GROUP=RESIN, NAME=matrix, MATERIAL=epoxy,TEMP=NO, MECH=HOM
Definition of HOMAT specific Commands¶
In the following section the HOMATspecific keywords are listed not alphabetically but in the order they are classically written by Mesh2HOMAT in the HOMAT command file PROBLEM_NAME.hin
.
Options introduced before the RVE description¶

/HEADING
The command file
PROBLEM_NAME.hin
always starts with the/HEADING
keyword. It has no specific parameters and permits the user to specify the intended homogenization simulation by providing a problem description. When generated by Mesh2HOMAT, the header section of the command file contents 6 lines, where the second line is left blank for the user to add a short problem description. The first line contains the date and time of the Mesh2HOMAT call. The 3^{th} until 6^{th} line outlines the total number of elements, nodes, new generated nodes and materials respectively.Note: It is important that the total number of nodes and new nodes indicating in the heading corresponds to the number of nodes specified in the
PROBLEM_NAME.hgeo
file. HOMAT checks this coincidence. Ex. 13 of the thermoelastic homogenization of a UD reinforced composite starts with following header:/HEADING Date: 20.11.19 Time: 12:45:24 Total number of Elements: 296 Total number of Nodes: 732 Total number of new Nodes: 76 Total number of Materials: 2

/PROBLEM
,TYPE =
type
,HOM_MOD =
hom_mod
,PER_MOD =
per_mod
This option defines the desired homogenization analysis and is the keyword with the largest impact on the predicted effective properties. The nature of the homogenization analysis is specified here via the parameter
TYPE
. The user has also to choose one of following three homogenization models: 3D (default), 2D plane strain or 2D plane stress model, via the parameterHOM_MOD
. The distinction between plane stress and strain is only relevant for thermoselastic homogenization. If a 2D thermal or fluid flow homogenization is desired,HOM_MOD = 1
is always set. The user must specify the nature of the imposed periodicity boundary conditions: 1D, 2D or 3D. Since HOMAT 6.0, 3D periodic BCs are assumed by default. TheTYPE
specifier can take the following values:TYPE=THERMAL
A purely thermal homogenization analysis is performed and only the effective thermal conductivities, capacities and densities are evaluated.TYPE=THERMEC
A full thermomechanical homogenization analysis is performed. In this case, a thermal homogenization and a thermoelastic one are realized successively in the same run. Effective thermal conductivities, capacities and densities are determined as well as the effective Hooke matrix, the effective thermal expansion coefficient and the generalized eigenstrain, if present.TYPE=MECHAN
A purely mechanical homogenization run is performed. In this case, only the Hooke matrix and its derivate effective properties (Young and shear moduli, Poisson coefficient,...) are evaluated but not the effective thermal expansion coefficient and the volumetric eigenstrain.TYPE=THELAS
A thermoelastic homogenization analysis is performed. In this case all effective thermoelastic properties (Hooke matrix and its derivate, effective thermal expansion coefficient, ..) are determined, but no thermal properties.TYPE=STOKES
A Stokes flow homogenization analysis is performed. The effective Darcy permeability tensor is evaluated, by applying unitary boost volume forces on the fluid domains.TYPE=LOCALI
A mechanical localization analysis is performed including a previous mechanical homogenization step in order to calculate the required, microscopic, periodic displacement fields x^{rs} on the RVE but the effective mechanical properties are not evaluated. The microstresses and strains are evaluated on the RVE for a given macrostrain state.TYPE=LOCATH
A thermoelastic localization analysis is performed including a necessary thermoelastic homogenization step. Total microstrains, including thermal expansion and volumetric eigenstrains, and microstresses are evaluated on the RVE for a given macrostrain state.
The
HOM_MOD
parameter can take the following values:HOM_MOD = 0
3D homogenization model. In this case, for the thermal analysis, three linear problems are solved. They correspond to the microscopic displacement fields \zeta^{X}, \zeta^{Y} and \zeta^{Z} respectively. In the case of mechanical homogenization, HOMAT establishes and solves six linear systems, corresponding to the microscopic displacement fields x^{rs} (x^{XX},\, x^{YY},\, x^{ZZ},\,x^{XY},\, x^{XZ},\, x^{YZ}) respectively. These fields are induced by the corresponding unitary macrostrain tensor \vec{E}_{rs}. In the case of flow homogenization, three special Stokes flow problems are solved providing the velocity fields \omega^{x}, \omega^{y} and \omega^{z} and the pressure fields p^{x}, p^{y}, p^{z} on the RVE respectively. The 3D homogenization model is specified by default.HOM_MOD = 1
2D plane strain homogenization model. In the case of thermal or flow homogenization, two linear problems are solved in the XY plane. They provide either the planar microscopic displacement fields \zeta^{X}, \zeta^{Y} in the thermal analysis or the planar velocity and pressure fields (\omega^{x}, \omega^{y} and \omega^{z} and p^{x}, p^{y}). 2D plane strain mechanical behaviour (\epsilon_{ZZ}=0 but \sigma_{ZZ} \neq 0, corresponding to a section through an infinite thick component) is defined in the mechanical homogenization leading to solve four linear problems, corresponding to the displacement fields x^{XX}, x^{YY}, x^{ZZ},x^{XY}.HOM_MOD = 2
2D plane stress homogenization model. In the case of thermal and flow homogenization, this model equals the previous one. No difference between plane strain and plane stress assumptions is made there. In the case of mechanical homogenization, the plane stress assumption (\sigma_{ZZ}=0 and \epsilon_{ZZ} =0, thin membrane behaviour) leads to solve only three linear systems, providing the displacement fields x^{XX}, x^{YY}, x^{XY}.
The
PER_MOD
parameter determines the nature of the applied periodic boundary conditions:PER_MOD = 1
1D periodic boundary condition is specified. In this case, two opposite free surfaces are linked. Classically, they are located normal to the X direction, but they can also be defined in other direction. The keyword/PERIODICITY
defines the pair of linked surfaces.PER_MOD = 2
2D periodic boundary conditions are defined. In this case two pairs of linked surfaces must be specified for an RVE of cubic or rhombic shape via the/PERIODICITY
option. Classically, the free surfaces in X direction are linked together as well as both free surfaces in Y direction. But, it is also possible to specify other linked pair of opposite surfaces, e.g. in X and Z directions or in Y and Z directions.
PER_MOD = 3
3D periodic boundary conditions are applied on the RVE. In this case, three pairs of linked surfaces must be specified via the/PERIODICITY
option.Note: At the moment, Mesh2HOMAT automatically generates the periodic boundary conditions only for RVEs of cubic or rhombic shape but not for unit cells of hexagonal or octahedral shape (see Fig. 10 and 11).
The
SOLVER
option specifies the adopted solver for the discretized linear global system.
SOLVER=ITER
An improved iterative BiCGstab solver is used, by default, to solve the linear systems of the mechanical, thermal or flow problems. This solver is coupled with a symmetric preconditioner and uses, since version 6.0, a sparse matrix storage format for the system matrix.Note: At the moment, the iterative BiCGstab solver is the only available solver for the thermal, mechanical and fluid flow homogenization problems.
In Ex. 16, the 2^{nd} part of the command file for the thermomechanical homogenization of a UD reinforced NiAl composite with Al_{2}O_{3} (sapphire) fibres is outlined. It is the part of the command file with all specific HOMAT keywords, which are introduced after the keywords describing the RVE model. The
/PROBLEM
keyword indicates there that successive thermal and thermoelastic homogenizations will be performed on a RVE with square shape on which 2D periodic boundary conditions are applied.
# /PROBLEM, TYPE=THERMEC, HOM_MOD=0, PER_MOD=2 # /DEF LAYER, TYPE=MULTI, NUMBER=1 LAYER=1, NAME=RVE, NB_DOM=ALL # /TEMP DOMAIN, NUMBER=2, NTEMP=2 DOMAIN=1, 800.0, 1100.0 DOMAIN=2, 800.0, 1100.0 # # Solver settings # /NUMERICS, ITMAX=2000, PRECISION=5.0D05, PENAL=1.0D+05, ITCRIT=1, OUTFREQ=20, EPSCOORD=1.0D02, TINY=1.0D16 # # Parameter definition for the homogenization module # /HOM PARA, MINTEMP=0.0, MAXTEMP=1400.0, CUBVAL=0.1, TREF=20.0 # /OUTPUT, FORMAT=VTK, TYPE=ASCII, IOP=2, DBG=YES, VOLAVRG=YES, VISLAV=YES T, U # /GAUSS INTEGRATION,* TYPE=HEX, N_X=2, N_Y=2, N_Z=2 # /GAUSS INTEGRATION*, TYPE=PENT, N_T=3, N_Z=2 /GAUSS INTEGRATION*, TYPE=TET, NUMBER=4 # /GAUSS INTEGRATION,* TYPE=QUAD, NUMBER=4 /GAUSS INTEGRATION,* TYPE=TRIA, N_T=3 #

/OUTPUT
This keyword permits the user to define the output of the considered homogenization or localization analysis. Sec. HOMAT Output details the different output possibilities for homogenization and localization analyses. The parameters of the keyword line are the following:

FORMAT
: this parameter is set by default toVTK
. In fact, the VTK format is used as standard exchange protocol between HOMAT and other programs, working on the RVE level, such as MICRESS or SphaeroSim. Moreover, HOMAT relies on the opensource, multiplatform visualization program Paraview to visualize the implicit fluxes/forces and periodic displacement, microstrain, microstress, velocity and pressure fields. 
TYPE
: the VTK output files can be eitherASCII
orBINARY
. TheBINARY
option allows to reduce the size of the result files and is the default setting if the HOMAT analysis is performed after a microstructure simulation with MICRESS or SphaeroSim. For other homogenization and localization runs, theASCII
option is adopted as the defaultTYPE
(see Ex. 16). 
DBG
: this parameter triggers additional control output either on the screen or in a log file, namedPROBLEM_NAME.hout
. By default, this parameter is set toYES
. The activation of this output option is recommended as it does not produce much overhead, but allows to follow the individual program steps in more detail. 
IOP
: this option controls the level of printed intermediate results, like the Hooke matrix per grain/material, the evaluated implicit forces, the element stiffness matrix and so on either on the screen or in thePROBLEM_NAME.hout
file. This output control parameter varies between 0 and 5 and a higher value produces more output. IfIOP=0
, the screen result output is strictly limited. This value corresponds to the default option for large scale applications. Note thatIOP=2
is a very frequently used choice and recommended if the user want to control the main intermediate output results. Values aboveIOP=2
are only needed to debug the program. 
ENG
: this parameter allows, ifENG=YES
, the direct output in the .mat file of the effective orthotropic Young and shear moduli and of the Poisson coefficients. This output is not classical ABAQUS one, but can be used by other macroscale simulation tools. By default,ENG=NO
. 
VOLAVRG
: ifVOLAVRG=YES
, then in the effective properties result file,PROBLEM_NAME.hres
, not only the effective thermal and/or thermoelastic properties are displayed; but also the corresponding mean volume averaged properties of the layers (monolayer homogenization) or of the whole RVE (multilayer homogenization). This way, the user can better assess the accuracy increase, induced by the correction term of the first order homogenization [cf. Eq. \eqref{eq:solu_3}]. Mesh2HOMAT writes in the template of the command filePROBLEM_NAME.hin
the default valueVOLAVRG=YES
. 
POLE
: this parameter allows to draw the pole figure in the plane XY of the lamella orientations perpendicular to the radial growth direction in each spherulite of semicrystalline thermoplastic microstructure. The pole figure are realized by stereographic projection of the orientations on the plane XY (see Fig. 12). 
VISLAV
this parameter permits the accurate visualization of forces and fluxes within Paraview. It allows to regenerate the original values of the forces and fluxes in the master nodes, which in the solved reduced system content not only its value but also the values of corresponding linked slave nodes. By default,VISLAV=YES
. 
LOC_AXIS
: ifLOC_AXIS=YES
, then the effective material properties are not only written in the structural RVE axis system but also in the specified local axis system in both result files:PROBLEM_NAME.hres
andPROBLEM_NAME.mat
. Of course, the local axis system has to be specified before via the/ORIENTATION
or the/MICROSIM
keyword. By default,LOC_AXIS=NO
. 
SHEAR_TENS
this parameter governs the output of the sheartension coefficients, which can be derived from a fully anisotropic effective Hooke matrix. IfSHEAR_TENS=YES
, these sheartension coefficients are calculated and printed in thePROBLEM_NAME.hres
file. In order to reduce the output of fully anisotropic Hooke matrices in this file,SHEAR_TENS=NO
by default.
In the data line of keyword
/OUTPUT
the user has to specify the nature of results fields on the RVE, which will be saved by HOMAT in the VTK format. A detailed description of the results files and their nomenclature is given in Sec. HOMAT Output. The following options are available:
ALL
: This option is selfexplanatory. All the evaluated microscopic fields are saved for possible visualization. By specifying this option, the user does not need to select the microscopic displacement and flux fields individually, results of the thermal homogenization via parameterT
or/and the displacement, microstrain, pressure and implicit force fields, results of the thermoelastic homogenization, via parametersU
,P
andE
. Both are saved in the case of a combined thermal and a thermoelastic homogenization. For a localization run, the microscopic strain and stress fields are both saved in separate files. For a Stokes flow homogenization, the boost forces, the velocity and pressure fields are saved in separated files. 
T
: In the case of thermal homogenization, the microscopic displacement and implicit flux fields are saved in separated files for visualization purpose. U
: In the case of thermoelastic or mechanical homogenization, the microscopic displacement fields and the implicit forces are saved in separated files for possible visualization.E
: In localization, mechanical as well as thermoelastic homogenization analyses, the variation of the microscopic engineering strain tensor over the considered layer or over the whole RVE are saved for visualization and/or posttreatment. In the case of localization, the equivalent strain is also saved for possible visualization.
S
: The components of the microscopic stress tensor and the equivalent stress is saved per layer or per RVE for visualization purpose in a localization run. 
P
: In presence of mixed finite elements (thermoelastic homogenization involving quasiincompressible materials or Stokes flow homogenization), the parameterP
indicates that the local pressure fields in these elements are saved for visualization purpose. 
V
: The components of the velocity fields, for Stokes flow problems, are saved per layer or RVE for visualization.
If the parameter
POLE=YES
is set, then the user has to specify via the additional parameterNB_SPHERU
on the third line (second data line) the number of spherulites for which the pole figured will be build. On the following line, the user specifies the list of spherulites for which a specific file containing the stereographic projection points of the lamella orientations, perpendicular to the radial growth direction, on the XY plane is built. The identification number of a spherulite is theNUMBER
value specified in the/SOLID SECTION
definition. In the spherulite example,mscal_PP
, we have introduced following output definition:/OUTPUT, FORMAT=VTK, TYPE=BINARY, DBG=YES, IOP=2, VOLAVRG=YES, ENG=NO, POLE=YES U NB_SPHERU=3 80, 103, 104
In Figure 13, the pole figure for the spherulite 103, shown in Fig. 12, is presented.


/UNITS
,GEOMETRY =
sfac
,LENGTH UNIT =
ME
CM
MM
MU
NM
,TEMP OFFSET =
offset
,ANGLE =
GRAD
RAD
This keyword permits the user to specify the desired unit system for the homogenization run. Following optional parameters can be specified on the keyword line:
GEOMETRY
This parameter allows the introduction of a scaling factor, which will be applied to all node coordinates. By default, no scaling is defined:GEOMETRY=1
.LENGTH UNIT
This parameter, linked to theGEOMETRY
parameter, can take either the valueME
for meter,CM
for centimeter,MM
for millimeter,MU
for micrometer andNM
for nanometer. It allows HOMAT to print, for example, the solid section volumes using correct units.MM
is the default value of this parameter.TEMP OFFSET
This parameter allows to change easily the temperature unit from grad Celsius to Kelvin. By default, the temperature is noted in °C. The offset is then: 273.15K.ANGLE
This parameter allows to specify the unit for angles: either as degree (ANGLE=GRAD
) or as radians (ANGLE=RAD
). By default, an angle is defined in degree.
In Ex. 13, following unit definitions are introduced:
/UNITS, GEOMETRY=1, TEMP OFFSET=273.15, ANGLE=GRAD
The geometry of the RVE is not scaled. Therefore, it was not necessary to specify the length unit as the default one is used. The angles are specified here in degree. If the user wants to use in the simulation microns as length unit and the original mesh is in millimeters, the following
/UNITS
parameters are specified:/UNITS, GEOMETRY=1.0E+3, LENGTH UNIT=MU, TEMP OFFSET=273.15, ANGLE=GRAD

/INCLUDE
,NAME =
name
The
/INCLUDE
keyword can be used at any time in the command file in order to include an external file that contains a part of the HOMAT command file, like thePROBLEM_NAME.hgeo
file, which describes the geometrical model of the heterogeneous material, whose effective properties will be evaluated. Note that, in order to reduce the size of the geometry file, Mesh2HOMAT writes, by default, a gzipped geometrical model. When a reference to an external file is encountered, HOMAT will unzip the file and process the data within the specified external file. At the end of this external file, HOMAT will automatically continue to interpret the calling file, which is generally the command filePROBLEM_NAME.hin
. The maximum number of nested/INCLUDE
operations is 4. This keyword has only one mandatory parameterNAME
. This parameter specifies the name of the external file. In Ex. 13, the RVE geometry and the thermophysical material properties of the fibre and matrix are specified externally and included via the option/INCLUDE
, as follows:/INCLUDE, NAME=ehz50_3D.hgeo.gz # /MATERIAL, NAME=VERRE, TYPE=ANISO /INCLUDE, NAME=verre_ani.hmat # /MATERIAL, NAME=EPOXY, TYPE=ANISO /INCLUDE, NAME=epoxy_ani.hmat

/MICROSIM
,NAME =
name
,TIME =
time
, [ORIE =
YES
NO
], [OUT_ORIE =
grain
]This keyword tells HOMAT that the microstructure to homogenize comes from a previous simulation with MICRESS. This option allows an easy and direct definition of the grain orientations of a polycrystal. In fact, the MICRESS problem name is known by HOMAT via the mandatory parameter
NAME
. The user has to provide in the working directory the ASCII fileMICRESS_NAME.TabO
containing the mean orientation per grain. As this text file contents orientation information for the whole microstructure evolution, it is essential to specify to HOMAT the MICRESS simulation time for which effective properties will be evaluated. This information is given via the mandatory parameterTIME
. Setting parameterORIE=YES
, indicates to HOMAT that, by default, the grains/regions of the microstructure receive their orientation directly via theTabO
file. Thus, the orientation definition via the keyword/ORIENTATION
and the parameterORIENTATION
of each solid section can be omitted. However, if this keyword is specified for a specific solid section, than its orientation definition overrides the default setting via theTabO
file.Remark: In the case of localization analysis with different structural axis system at micro and macroscale, the
/MICROSIM
keyword must be specified before the keyword/ORIENTATION
in the.hin
command file, if the user wants to specify the grain orientations easily via theTabO
file.In the presence of a 2D MICRESS microstructure, the parameter
ROTX=90
must be specified. It defines the rotation angle between the MICRESS axis system \vec{E}_{X}, \vec{E}_{Z}, \vec{E}_{Y} to the HOMAT structural axis system e_{X}, e_{Y}, e_{Z} (see Fig. 8). In the 2D case, HOMAT rotates, first, each vector or tensor quantity from the local crystal axis system to the MICRESS axis one and then from this system to its own reference system. This parameter has no relevance in 3D.The new optional parameter
OUT_ORIE
permits to specify the name of the grain orientation which specifies the local axis system for the output of the effective material properties in this local axis system. Only one grain can be selected for this additional output and a definition via the/ORIENTATION
keyword is obsolete.For the austenite to ferrite phase transformation of FeCMnSi steel, discretized in 3D and detailed further in the example manual, the user introduces the following parameters:
/MICROSIM, NAME=dim1.5, TIME=42.0, ORIE=YES, OUT_ORIE=grain2
The MICRESS name corresponds to the driving file name. Here, the microstructure after 42sec, containing 49.46% ferrite grains, is homogenized and all grain orientations are directly provided via the TabO file. Moreover, the effective thermoelastic properties are additionally written in the local axis system of
grain2
.
Options introduced after the RVE model description¶

/DEF LAYER
,TYPE =
MONO
MULTI
,NUMBER =
number
,LAYER =
,NAME =
name
,NB_DOM =
dom
In presence of a RVE build of different material layers, like the cooled multilayer build of a substrate of superalloy CMSX4, a MCrAlY bondcoat and a yttrium stabilized thermal barrier coating (see Fig. 14), the keyword
/DEF LAYER
allows to select either multilayer or monolayer homogenization via the keyword parameterTYPE
. Multilayer homogenization predicts the effective properties of the whole RVE; whereas in the case of monolayer homogenization, the user can specify individually the layers, which are homogenized in an analysis. Of course, in presence of only one material layer, both homogenization types are identical. Note that the multi or monolayer homogenization option can be specified in combination with any kind of homogenization analysis.
TYPE=MULTI
: multilayer homogenization
By default,TYPE
is set toMULTI
. In this case, the other keyword parameters are easily predefined by Mesh2HOMAT, like shown as following:/DEF LAYER, TYPE=MULTI, NUMBER=1, LAYER=1, NAME=RVE, NB_DOM=ALL
Only one layer is defined (
NUMBER=1
or0
), the number of the consideredLAYER
is always 1 and the layerNAME
is set toRVE
and the parameterNB_DOM
, specifying the number of considered domains or grains is set toALL
. 
TYPE=MONO
: monolayer homogenization
The number of layers considered by the homogenization run is introduced via parameterNUMBER
. If you desire the effective monolayer properties of each layer, the user introduces the list of the concerned material domains, as shown below for the cooled multilayer plate of Fig. 14:/DEF LAYER, TYPE=MONO, NUMBER=3 LAYER=1, NAME=TBC_HO, NB_DOM=6 1, 2, 3, 4, 5, 16 LAYER=2, NAME=BOND_HO, NB_DOM=6 6, 7, 8, 9, 10, 17 LAYER=3, NAME=SUBSTRATE, NB_DOM=6 11, 12, 13, 14, 15, 18
The number of materials domains per layer is specified via the parameter
NB_DOM
. Each layer of the cooled plate has 6 domains, the solid and 5 cooling holes. On the next line, the user has to specify, via their identification number, the list of the domains which belong to the layer. The identification number of a material domain is theNUMBER
value specified in the/SOLID SECTION
definition.
Remark: At the moment, HOMAT v6.0 allows to perform monolayer homogenization on layers of a multilayer RVE with 2D or 1D periodic boundary conditions. 3D periodic BCs on individual layers are not allowed.


/TEMP DOMAIN
,NUMBER =
number
, [NTEMP =
ntemp
]This keyword allows the definition of the domain temperature: either this temperature can vary inside each domain or the temperature is constant per domain. In the 1^{st} case a temperature field is provided, coming from a previous analysis via an input file (see option
/INPUT FILES
). In the 2^{nd} case, the user has to specify, for each domain, the temperatures at which homogenization will be performed in one run. This keyword requires the definition of the parameterNUMBER
, which specifies the nature of the temperature field.
In the first case, the user specifies only following keyword line:
/TEMP DOMAIN, NUMBER=0
Parameter
NUMBER
equals zero is directly interpreted by HOMAT that all domains have no constant temperature and the required local temperature field is imported via an input file. Only one homogenization simulation is performed with this set temperature profile on the RVE. 
In the second case, the user has, at first, to define via parameter
NUMBER
the number of material domains or grains affected by the present temperature definition. Per domain, a special entry line has to be introduced, as shown below. Further, the user must introduce the number of different temperatures, specified for each concerned material domain or grain, via the keyword parameterNTEMP
. IfNTEMP
differs from one, HOMAT will perform automaticallyNTEMP
homogenizations in one simulation run. This option is very interesting in applications where the temperature dependence of the effective properties is investigated. This feature allows to performNTEMP
homogenizations in a single run. At the moment, HOMAT limitsNTEMP
to 9 simultaneous homogenizations. An example is given below:/TEMP DOMAIN, NUMBER=2, NTEMP=3 DOMAIN=1, 27., 100., 200. DOMAIN=4, 27., 120., 220.
The keyword line states that two domains with three different temperature sets are specified. Per domain, the three set temperatures are introduced on each data line: The domains are identified by the parameter
NUMBER
of the/SOLID SECTION
keyword. Of course, the temperatures can be different, per domain/grain. HOMAT will perform a first homogenization at 27°C for both domains, then a 2^{nd} homogenization with domain 1 at 100°C and domain 4 at 120°C. Finally, it performs a third homogenization with domain 1 at 200°C and domain 4 at 220°C. This simple example shows the great flexibility of defining uniform temperatures per domain for different homogenization analyses. In order to condense the required input, mainly in presence of a polycrystal with several grains set to constant temperatures, the parameterNUMBER
can be set toALL
. In this case, the definition for each domain separately can be replaced by one common set of temperatures:/TEMP DOMAIN, NUMBER=ALL, NTEMP=2 27., 100.
In this example, two homogenization runs are specified with all domains set to 27.0°C and 100.0°C respectively. For a MICRESS microstructure, Mesh2HOMAT sets by default,
NUMBER=ALL
. 

/INPUT FILES
,INTEMP =
intemp
MACRO =
macro
, [ORIENTATION =
orientation
]This keyword permits to specify the name of different input files coming from linked, previous simulations, like a multiphase field analysis of metallic alloys with MICRESS. At the moment, only two parameters can be specified:
INTEMP
andMACRO
. They are briefly described below:/INPUT FILES, INTEMP=ReX_orien.vtk
The user introduces with this command line the file name of a temperature distribution over the RVE as initial condition for the homogenization run. This set temperature field may vary locally inside each domain/grain. The temperatures are there nodal values, specified via the POINT data format (see VTK documentation).
/INPUT FILES, MACRO=/your/path/to/U+O_forming.txt
The parameter
MACRO
defines the path to the ASCII file containing the macrostrain results of the linked thermomechanical macrosimulation. A python script, namedaba2hom
, extracts these macro strains from the ABAQUS results database (OBD file ^{6}) in critical areas of the analyzed structure. HOMAT uses these strain states to perform localization runs. The different contents of this file are outlined in Sec. Transfer of macroscale results from ABAQUS to HOMAT, dedicated specially to this topic.Note: The user can introduce either the local file name with extension in the working directory or the global name including the total path name, as outlined above. If the structural axis system in HOMAT differs from the used structural axis in the macroscale analysis, the user must specify the optional parameter
ORIENTATION
on the keyword line:/INPUT FILES, MACRO=Critical_orth_lin.txt, ORIENTATION=rot_macro
In this case, the orientation definition
rot_macro
allows HOMAT to rotate the macrostrain results of ABAQUS in the structural axis system of the microstructure in order to perform the localization analysis.

/NUMERICS
This keyword specifies numerical parameters for the iterative BiCGStab solver ^{17}. Mesh2HOMAT has endowed these parameters with a default value, which the user may adapt to his specific problem via the present command. Following parameters can be specified via this keyword:

ITMAX
The number of maximum allowed iterations for the iterative BiCGStab solver to converge to the unknown periodic displacement fields \zeta^{r} or x^{rs} or velocity fields \omega^{r}. If the number of maximum iterations is reached without fulfilling the convergence criterion, HOMAT takes the displacement field corresponding to the smallest residuum as the result. By increasing the maximum iteration number, the solver accuracy is increased but, in the same time, the CPU time increases. By default,ITMAX
is set to2000
. 
PRECISION
This parameter specifies the value, which the convergence criterion of the iterative BiCGStab must undershoot in order to converge. By default, in conjunction with the relative convergence criterion (ITCRIT=1
) this convergence parameter is set toPRECISION=5.10^{5}
. 
ITCRIT
This option specifies the adopted convergence criterion of the BiCGstab solver. By default, the relative convergence criterionITCRIT=1
is specified\begin{equation} \label{eq:relconv} \mathbf{\text{crit}}\left( \mathbf{x}^{i} \right) = \frac{\left \mathbf{A}\mathbf{x}^{i}  \mathbf{b} \right^{2}}{\left \mathbf{A}\mathbf{x}^{0}  \mathbf{b} \right^{2}} < \mathtt{PRECISION}\,, \end{equation}where \vec{A} is the symmetric system matrix and \vec{b} the right hand side vector of the linear system \vec{A}\cdot\vec{x} = \vec{b} and x^{0} and x^{i} are the initial guess of the displacement field and its expression at iteration i respectively. If
ITCRIT=2
, following absolute convergence criterion is adopted\begin{equation} \label{eq:absconv} \mathbf{\text{crit}}\left( \mathbf{x}^{i} \right) = \ \left \mathbf{A}\mathbf{x}^{i}  \mathbf{b} \right^{2} < \mathtt{PRECISION} \,. \end{equation} 
OUTFREQ
This option defines the frequency of printing intermediate solver results. The default valueOUTFREQ=20
means that after 20 iterations, HOMAT prints the value of the actual convergence criterion, the norm of the actual residuum: \vec{r}^{i} = \vec{A}\cdot\vec{x}^{i}  \vec{b}, the mean displacement or velocity and the best iteration. EPSCOORD
This numerical parameter is used in the implementation of periodic BCs to determine if a node has a direct correspondence on the master surface. This numerical parameter specifies a relative precision, which is multiplied by the mean element size of the slave face in order to define the local absolute precision. HOMAT ensures that the absolute precision is always smaller than the smallest element size on the considered free surface. By default,EPSCOORD=1.D02
. 
IPRE
This parameter specifies the adopted preconditioning scheme.IPRE=1
means that a diagonal matrix \vec{D} built of the square roots of the diagonal terms of system matrix \vec{A} is applied in symmetric way to the system matrix \vec{A}:\begin{equation} \label{eq:pre1} \vec{D}\cdot\vec{A}\cdot\vec{D}\cdot\vec{x} = \vec{D}\cdot\vec{b}\cdot\vec{D}\,,\text{ with } D_{ii} = \frac{1}{\sqrt{\leftA_{ii}\right}} \end{equation}In the case
IPRE=2
the precondition matrix \vec{D} is built from the absolute value of each diagonal term of matrix \vec{A} and only a lefthand side preconditioning is performed:\begin{equation} \label{eq:pre2} \vec{D}\cdot\vec{A}\cdot\vec{x} = \vec{D}\cdot\vec{b}\,,\text{ with } D_{ii} = \frac{1}{\leftA_{ii}\right} \end{equation}Numerical investigation has shown that the new symmetric preconditioning outperforms the leftside one, used in HOMAT v5.x versions. Therefore,
IPRE=1
is the default option. 
TINY
This parameter defines a small number, which is used throughout HOMAT as numerical equal zero. By default,TINY=1.0E16
. 
ADMIS
This parameter permits to specify the max. admissible increase of the mean solution (displacement or velocity) between two successive iterations: If the actual mean solution exceeds its admissible increase, a warning message is written in the log filePROBLEM_NAME.hout
. Mainly in the first iterations such a warning appears. If such message appears more than five times during the solution of a specific linear system, a detailed analyzed of the problem and its definition is recommended. By default,ADMIS=200
.


/HOM PARA
This keyword allows the specification of different global parameters, which are used throughout the homogenization analyses.

MINTEMP
defines the lowest temperature (in °C) used in the definition of all temperature dependent thermophysical properties. By default,MINTEMP=0
. 
MAXTEMP
defines the maximum temperature (in °C) used in the definition of all temperature dependent thermophysical properties. By default,MAXTEMP=1400
. 
CUBVAL
defines the precision value used in the criterion for detecting cubic symmetry in the RVE coordinate system. If\begin{equation*} \frac{\left H_{\text{ii}}  H_{\text{mean}} \right}{H_{\text{mean}}} < \mathtt{CUBVAL}\,, \end{equation*}the equivalent material is considered cubic^{27}. Here H_{\text{mean}} = \frac{H_{\text{ii}}}{3}, i = 1,3 is the mean value of the diagonal Hooke matrix terms H_{ii}. In this case, the effective Young's modulus E, the elastic constants C_{11}, C_{12} and C_{44} and the anisotropic Zener coefficient A^{Z} [see Eq. \eqref{eq:eff_zener}] are provided in addition to its effective orthotropic Young's and shear moduli E_{i}, G_{ij} and Poisson coefficients \nu_{ij} in the
PROBLEM_NAME.hres
output file. If\begin{equation*} \frac{\left H_{ii}  H_{\text{mean}} \right}{H_{\text{mean}}} > \mathtt{CUBVAL}\,, \end{equation*}the equivalent material is orthotropic. Only its effective Young moduli, E_{i}, per structural direction i, the shear moduli G_{ij} per plane ij and the Poisson coefficients \nu_{ij} are given in the
.hres
result file. By default,CUBVAL
is set to0.01
. 
HTR_IDEAL
specifies the value of an ideal heat transfer between two materials. By default, this value is set to0.01
W/(mm^{2}K). This value is used for all interfaces with no specified ideal heat transfer coefficient. 
TREF
This parameter specifies the reference temperature in °C for the homogenization and localization runs. This reference temperature plays an important role in the definition of eigenstrains via molar volume of the different phases/materials, as it defines the initial state without any eigenstrain in the heterogeneous material. In the thermoelastic localization analysis this temperature specifies also the reference strainless state. By default,TREF=20
. 
FOAMBOX
IfFOAMBOX=YES
, then the user has the possibility to define a box inside the RVE where the velocities inside the porous medium is used to express the Darcy permeability. Velocity values outside this box are ignored in the permeability evaluation. This way, the impact of the periodic boundary condition on a unit cell of a material with random microstructure and of delicate pore cutting, like for opencell metallic foams, is reduced in the permeability prediction ^{11}. By default, this parameter is set toNO
. IfFOAMBOX=YES
, the user has to specify on the next line the extension of the foam box in the adopted unit system:/HOM PARA; MINTEMP=0, MAXTEMP=1400., FOAMBOX=YES 0.1, 0.2, 0.2

RANDOM
This parameter specifies the seed used by the random number generator. Random numbers are defined, for example, by the spherulite model to generate, for each spherulite, the initial orientation of the bilamella of crystalline and amorphous phases in the plane perpendicular to the radial growing direction. If no random seed is specified by the user, HOMAT sets it to256
.


/GAUSS INTEGRATION
The finite element mesh of the RVE is composed of hexahedral, pentahedral and/or tetrahedral elements. The variation of the field (temperature, ..) inside an element is approximated either via a linear or quadratic shape function. The evaluation of quantities over the element volume, like the element stiffness matrix, is performed numerically via the GaussLegendre method, which requires, per element type, the definition of integration points, named Gauss points. HOMAT offers the possibility to choose between a full or a reduced integration schemes ^{10}. For each element type, the default definition of the Gauss integration scheme set by Mesh2HOMAT is as follows:
 linear and quadratic hexahedron: 2 \times 2 \times 2 Gauss points, which means that 2 Gauss points are defined per structural direction, as shown in Fig. 15(a).
 linear and quadratic pentahedron: 3 \times 2, means 3 Gauss points in planes parallel to both triangular faces and 2 Gauss points in the element thickness, defining the Z position of these two planes [see Fig. 15(b)].
 linear and quadratic tetrahedron: 4 Gauss points, located as shown in Fig. 15©.
 linear and quadratic quadrangular element: 2 \times 2 Gauss points, located as shown in Fig. 15(d).
 linear and quadratic triangular element: 3 Gauss points, located as shown in figure Fig. 15(e).
As reduced integration can lead to the development of hourglass modes (known as zero energy modes) ^{10} and no hourglass control is implemented in HOMAT, reduced integration, like using only one Gauss point for a linear hexahedral element, is extremely delicate to adopt. Therefore, the modification of the number of default integration points is reserved to advanced users, having profound knowledge about the finite element method. Per element type, the user has following possibilities for integration:
 linear hexahedron element: 1 \times 1 \times 1 or 2 \times 2 \times 2 Gauss points;
 quadratic hexahedron element: 2 \times 2 \times 2 or 3 \times 3 \times 3 Gauss points;
 linear pentahedron element: N_{T} = 1, 3 or 4 and N_{Z} = 1 or 2, where N_{T} are the number of GP defined per triangular face and N_{Z} the number in the thickness direction;
 quadratic pentahedron element: N_{T} = 3, 4, 6 or 7 and N_{Z} = 2 or 3 Gauss points;
 linear tetrahedral element: 1 or 4 Gauss points;
 quadratic tetrahedral element: 4 or 5 Gauss points;
 linear quadrangular element: 1 \times 1 or 2 \times 2 Gauss points
 quadratic quadrangular element: 2 \times 2 or 3 \times 3 Gauss points
 linear triangular element: 1 or 3 Gauss points
 quadratic triangular element: 3 or 4 Gauss points
Note: Independent of the element type, the new stabilized mixed elements (equal order, LiHe, mini or bubble) have the same Gauss point definition as the corresponding linear element version.
To specify the integration scheme in HOMAT, the user must define on the keyword line the concerned element type:
TYPE=HEX
,TYPE=PENT
orTYPE=TET
for a hexahedral, pentahedral or tetrahedral volume element, respectively, orTYPE=QUAD
orTRIA
for a quadrangular or triangular surface element. Then, he has the choice to specify the number of integration points directly via the parameterNUMBER
or by specifying the number of integration points per direction. For example, he can specify: N_{X}=2, N_{Y}=2 and N_{Z}=2. Note that the intrinsic directions of each element type are defined by the element node numbering (see Figs. 5 and 6). Examples 17 and 18 illustrate different equivalent Gauss point definitions for linear elements:# /GAUSS INTEGRATION, TYPE=HEX, N_X=2, N_Y=2, N_Z=2 # /GAUSS INTEGRATION, TYPE=PENT, N_T=3, N_Z=2 /GAUSS INTEGRATION, TYPE=TET, NUMBER=4 # /GAUSS INTEGRATION, TYPE=QUAD, NUMBER=4 /GAUSS INTEGRATION, TYPE=TRIA, N_T=3 #
# /GAUSS INTEGRATION, TYPE=HEX, NUMBER=8 # /GAUSS INTEGRATION, TYPE=PENT, NUMBER=6 /GAUSS INTEGRATION, TYPE=TET, N_T=4 # /GAUSS INTEGRATION, TYPE=QUAD, N_X=2, N_Y=2 /GAUSS INTEGRATION, TYPE=TRIA, NUMBER=3 #

/EIGENSTRAIN
,NUMBER =
number
This keyword must be specified in thermoelastic homogenization, if the user wants to take isotropic volumetric eigenstrains into account. Eigenstrains can be introduced in HOMAT in two different ways:

In presence of a temperature dependent molar volume of each involved material or phase, the user has only to specify on the keyword line via the parameter
MATRIX
, the material or phase name of the matrix. All other materials or phases are then considered as inclusions, producing a volumetric eigenstrain. The reference temperature,TREF
, specified by the/HOM PARA
keyword, defines the temperature T_{\text{ref}} at which it is assumed that no eigenstrains exist in the heterogeneous material. Thus, it defines the starting temperature for eigenstrain development in all inclusions via expression:\begin{equation} \label{eq:eigenstrain} \epsilon_{\text{eigen}} = \frac{\sqrt[3]{M^{\text{vol}}_{\text{incl}}(T)}  \sqrt[3]{M^{\text{vol}}_{\text{matrix}}(T)}}{\sqrt[3]{M^{\text{vol}}_{\text{matrix}}(T)}} + \frac{\sqrt[3]{M^{\text{vol}}_{\text{matrix}}(T)}  \sqrt[3]{M^{\text{vol}}_{\text{matrix}}(T_{\text{ref}})}}{\sqrt[3]{M^{\text{vol}}_{\text{matrix}}(T_{\text{ref}})}} \end{equation}This easy definition is recommended to be used in combination with a MICRESS simulation, for example, by solid state transformation of a steel microstructure, as the user has only to specify the matrix material.  constant eigenstrain can be defined for each material or phase by introducing the number of inclusion materials/phases presenting a constant eigenstrain via parameter
NUMBER
on the keyword line. Then on the next lines, the user has to specify, per inclusion, the name of the corresponding phase via the parameterPHASE
and the constant value of the isotropic eigenstrain. Ex. 19 illustrates the definition of constant eigenstrain per inclusion for a steel alloy where the austenite phase is the matrix and ferrite and pearlite are considered as inclusions.# /EIGENSTRAIN, NUMBER=2 PHASE=FERRITE, 4.3135E03 PHASE=PEARLITE, 3.2137E03 #
Note: This definition neglects, at the moment, the anisotropic character of the eigenstrain of pearlite composed of lamellae of ferrite and cementite.


/LOCA
,NUMBER =
number
This option defines the number of critical locations in the sample at the macroscale, for which a localization run will be performed on the RVE. On the keyword line the user has to specify the number of desired localization analysis via parameter
NUMBER
. Then, on the next lines and, for each localization analysis, the user specifies the element number of the macromodel via the parameterELE
and the Gauss point number via parameterGP
. IfGP=0
, then all Gauss points of the element are considered. Then, depending of the element type at the macroscale, defined in the transfer file (described in Sec. Transfer of macroscale results from ABAQUS to HOMAT),NGP
localization analyses will be performed. Moreover, if the parameterELE=0
, then all elements of the macrosimulationELSET
, defined in the transfer file via command (/INPUT FILES, MACRO=NAME_MACROSIM.txt
) are selected. For example,ELE=0, GP=1
, means that a localization run is performed for all elements of the transferred macro ELSET with the strain state taken always from the first Gauss point of each element.Ex. 20 illustrates the definition of 6 successive localization runs after Oforming of a pipeline tube (plane strain finite element analysis), which will be performed in the same HOMAT run.
# /LOCA, NUMBER=3 ELE=1297, GP=1 ELE=1253, GP=0 ELE=1248, GP=3 #
On the keyword line, parameter
NUMBER
gives the number of lines with specific localization input. For macroelements1297
(GP=1
) and1248
(GP=3
) a localization run will be performed as well as for all Gauss points of macroelement1253
. As a plane strain forming simulation has been achieved at the macroscale with linear quadrangular elements having 4 Gauss points, the definition of Ex. 20 will result in 1+4+1=6 localization runs. If in the macrosimulation anELSET=critical_elements
is defined, which involves the elements1297
,1253
,1248
, then the following input leads to three localization runs: for each element of theELSET
the strain state in the first Gauss point is selected to be applied as boundary condition on the RVE of the steel microstructure:/LOCA, NUMBER=1 ELE=0, GP=1

/SPHERU CENTER
,ANGSPH =
angsph
,NUMBER =
number
If the spherulite model (
TYPE=SPHERU
) is selected in the keyword/MATERIAL
for several material domains, the user must then also specify the spherulite centres of the semicrystalline microstructure. The microstructure evolution program SphaeroSim provides these centres to Mesh2HOMAT, which put them in a specific file, namedPROBLEM_NAME.hsph
, which can be read by HOMAT. On the keyword line, the mandatory parameterNUMBER
has to be specified. It defines the total number of spherulite centres, specified in the RVE for the considered homogenization analysis. On the following lines and per spherulite, the user enters its identification number (ID), followed by the three coordinates of its centre. Note that the spherulite ID must corresponds to the number introduced in the/SOLID SECTION
definition. A second, facultative parameterANGSPH
on the keyword line defines the rotation angle increment (in grad) between two successive axis locations in the plane perpendicular to the radial growth direction. This angle increment allows for a spherical spherulite to have an effective quasiisotropic elastic behaviour. Numerical experiments on a single spherical spherulite have shown that this angle lies between 60 deg. and 90 deg. By default,ANGSPH=80
deg. More details about the definition of this angle increment can be found in ^{3}./SPHERU CENTER, ANGSPH=70, NUMBER=756 1, 0.017, 0.009, 0.001 ... 736, 0.065,0.043, 0.099
In this example, a RVE of extension 0.1 mm^{3} with 756 solidified spherulites is analysed. In the spherulite model an angle increment of 70° deg is adopted here. The coordinates of the spherulite centres are specified here in mm.
HOMAT Output¶
In the following section the different HOMAT result files and their nomenclature are introduced. Sec. Effective properties results concerns the description of the files containing the predicted effective properties of the considered heterogeneous material. In Sec. Field Results the nature and nomenclature of the VTK results files of a homogenization analysis are described in detail. These files contain field variables, which vary in the RVE of the heterogeneous material. In Sec. Localization Results the nature and nomenclature of the field output files of a localization analysis are outlined. Eventually, in Sec. Log File Output a brief overview of the principal control outputs in the log file, named either PROBLEM_NAME.hout
for homogenization or PROBLEM_NAME.lout
for localization is given.
Effective properties results¶
In the case of a monolayer or multilayer homogenization analysis, HOMAT provides always two results files: the first one is named PROBLEM_NAME.hres
. It is an ASCII file containing the evaluated effective properties, presented in a documented and humanreadable way. The second file, named PROBLEM_NAME.mat
, contains the effective properties, written in the ABAQUS format for subsequent structural analysis on the macroscale.
Thermal homogenization results¶
For each selected layer or for the whole multilayer structure, HOMAT provides the calculated effective, orthotropic thermal conductivity tensor, expressed in the structural axis system and usually in W/(mm K). Moreover, the effective heat capacity [in J/(mm^{3} K)], density [in g/(mm^{3})] and specific heat [in J/(g K)], obtained by volume averaging (or homogenization of order zero), are printed in the result file .hres
. Additionally, if desired by the user, HOMAT writes also the effective thermal conductivity tensor in the adopted local axis system. Ex. 22 shows the thermal homogenization result file .hres
for a layer of a unidirectional (UD) composite with rectangular fibre shape, defined by Dvorak as a benchmark test ^{18} for his homogenization method. The rectangular fibre has a conductivity of 10^{2} W/(mm K); while the matrix has a larger conductivity of 1.0 W/(mm K). With his numerical approach, Dvorak predicts the following orthotropic effective thermal conductivities in the crosssection: \lambda_{1} = 0.671 W/(mm K) and \lambda_{2} = 0.425 W/(mm K). HOMAT provides with a coarse mesh of 12 \times 12 \times 2 linear hexahedral elements the effective properties given in Ex. 22. The result file .hres
starts always by giving the name of the present homogenization analysis, here hexas_scale
, followed by the number of specified layers (here one) and the number of homogenizations, here also one. Then, the effective orthotropic thermal conductivity tensor is written in the structural axis system.
 
 HOMOGENISATION RESULTS  EFFECTIVE PROPERTIES 
 
 SIMULATION: hexas_scale_1 
 
 number of layers : 1 
 number of homogenisations : 1 
 

  
 heat conductivity  k_i1  k_i2  k_i3 
    
 k_1j  6.776877E01  1.803914E12  
 k_2j  1.109493E11  4.272344E01  
    
 mean value  Dk_x  Dk_y  Dk_z 
    
 5.524610E01  22.67  22.67  0.00 

    
 heat capacity  rhocp   
    
  6.500015E02   
    
 Density  dens   
    
  2.062501E03   
    
 specific heat  cp   
    
  2.750264E+01   
    


 heat cond. vol. aver. kav_i1  kav_i2  kav_i3 
    
 kav_1j  7.502492E01  0.000000E+00  
 kav_2j  0.000000E+00  7.502492E01  
Due to the rectangular shape of the fibres, the equivalent layer presents a pronounced orthotropic thermal behaviour with \lambda_{1} = 0.67769 W/(mm K) and \lambda_{2} = 0.42723 W/(mm K). In spite of the coarse mesh, these values differ from Dvorak's predictions only by about 1% for λ_{1} and by 0.62% for \lambda_{2}. The next line provides the mean value and the relative deviation for each component from this mean value. A deviation of ± 22.67 % from the mean, isotropic conductivity is computed, which emphasizes the significant orthotropic character of the equivalent material. Effective heat capacity, density and specific heat are also provided. As the output parameter VOLAVRG=YES
is set, the volume averaged thermal conductivity is written in the output file. As the volume averaging procedure takes only the volume fraction of the fibre and not its rectangular shape into account, it predicts an erroneous isotropic equivalent behaviour! This simple example shows the advantage of a first order asymptotic homogenization scheme in comparison to the volume averaging method. Ex. 23 shows the contents of the .mat
file for a further thermal analysis at the macroscale with ABAQUS. As no name for the investigated composite is provided, HOMAT denotes the equivalent layer simply as a RVE material. The user can change this default name, by editing the .mat
file.
,\*,\*
,\*,\*
,\*,\* NUMBER OF SPECIFIED MATERIALS: 1
,\*,\*
,\*MATERIAL, name = RVE
,\*,\*
,\*DENSITY
2.062501E03
,\*CONDUCTIVITY, TYPE=ORTHO
6.776877E01, 4.272344E01, 5.524610E01
,\*,\*
,\*SPECIFIC HEAT
2.750264E+01
Thermal and thermoelastic homogenization results¶
In order to illustrate the versatility of HOMAT, we predict in the same homogenization run effective thermal and mechanical properties of a metal matrix composite at different selected temperatures. The example of a FG75 NiAl matrix reinforced by UD sapphire fibres of cylindrical shape is considered here. It is assumed that the fibres have a volume content of 30%. The ASCII output result file .hres
starts, after the definition of the problem name, with the print of the thermal homogenization results, followed by the thermoelastic ones. Example 24 shows the thermal part of the .hres
result file.
 
 HOMOGENISATION RESULTS  EFFECTIVE PROPERTIES 
 
 SIMULATION: fas30_thmec_ani 
 
 SIMPLE PROBLEM DEFINITION 
 
 number of layers : 1 
 number of homogenisations : 2 
 

  
 heat conductivity  k_i1  k_i2  k_i3 
    
 800.00 
 k_1j  3.434308E02  6.504486E10  0.000000E+00 
 k_2j  1.710552E11  3.434307E02  0.000000E+00 
 k_3j  6.310036E11  4.336384E11  3.937139E02 
    
 1100.00 
 k_1j  3.133230E02  8.107663E10  0.000000E+00 
 k_2j  5.581284E11  3.133228E02  0.000000E+00 
 k_3j  7.202623E11  5.052411E11  3.769298E02 
    
    
 mean value  Dk_x  Dk_y  Dk_z 
    
 3.601918E02  4.65  4.65  9.31 
 3.345252E02  6.34  6.34  12.68 

    
 heat capacity  rhocp   
    
 800.00  3.076047E03   
 1100.00  3.532198E03   
    
 Density  dens   
    
 800.00  4.742357E03   
 1100.00  4.727331E03   
    
 specific heat  cp   
    
 800.00  6.559209E01   
 1100.00  7.480665E01   
    
For each thermophysical quantity, effective orthotropic or isotropic properties are denoted at both defined temperatures, here at 800°C and 1100°C respectively. The cylindrical shape (see Fig. 4) induces an orthotropic effective behaviour. Indeed, deviations from the mean thermal conductivity of 4.65% (\lambda_{11} and \lambda_{22}) and of 9.31% for \lambda_{33} in fiber direction are predicted respectively at 800°C. Slightly higher deviations are predicted at 1100°C. The mechanical result part (see Ex. 25) starts by providing the calculated effective Hooke matrix of the UD reinforced composite, before its symmetrization. This way, the user can easily check how far the calculated Hooke matrix is from the expected symmetric one and obtains an error measure for the impact of the numerical discretization. As the sapphire fiber is monoclinic with only one symmetry plane perpendicular to the fibre, the effective Hooke matrix of the UD reinforced composite is fully anisotropic with extensionshear coupling terms. In order to derive effective engineering properties, the calculated Hooke matrix is symmetrized and inverted. The deduced effective shear and Young moduli as well as the Poisson coefficients of the composite layer are given in Example 26.
 Hooke matrix 6*6  H_i1  H_i2  H_i3 
 800.00 
 H_1j  3.001139E+05  1.328843E+05  1.229083E+05 
 H_2j  1.328979E+05  3.000802E+05  1.228837E+05 
 H_3j  1.195747E+05  1.195287E+05  3.190173E+05 
 H_4j  3.502496E03  2.038315E03  1.421406E04 
 H_5j  2.420831E04  2.026955E04  2.636610E08 
 H_6j  2.955201E+03  2.814527E+03  4.557293E02 
    
 part II  H_i4  H_i5  H_i6 
    
 H_1j  3.836078E03  7.183098E04  2.301368E+03 
 H_2j  1.579003E03  3.476790E04  2.301262E+03 
 H_3j  1.188207E04  1.084709E05  1.058137E01 
 H_4j  8.040503E+04  1.258861E+03  5.246791E04 
 H_5j  2.525469E+03  8.446847E+04  5.621908E04 
 H_6j  1.232085E04  3.716822E04  8.550840E+04 
 Shear moduli  G_xy  G_xz  G_yz 
    
 800.00  8.036740E+04  8.442893E+04  8.542898E+04 
 1100.00  7.131272E+04  7.869988E+04  7.938645E+04 
    
 mean value  DG_xy  DG_xz  DG_yz 
    
 8.340844E+04  3.65  1.22  2.42 
 7.646635E+04  6.74  2.92  3.82 

    
 Young moduli  E_x  E_y  E_z 
    
 800.00  2.242213E+05  2.242187E+05  2.511521E+05 
 1100.00  2.018505E+05  2.018355E+05  2.332937E+05 
    
 mean value  DE_x  DE_y  DE_z 
    
 2.331973E+05  3.85  3.85  7.70 
 2.123265E+05  4.93  4.94  9.87 

    
 Poisson coefficients  nu_xy  nu_xz  nu_yz 
    
 800.00  3.422293E01  2.800512E01  2.798895E01 
 1100.00  3.459106E01  2.791244E01  2.791171E01 
    
 mean value  Dnu_xy  Dnu_xz  Dnu_yz 
    
 3.007234E01  13.80  6.87  6.93 
 3.013840E01  14.77  7.39  7.39 

    
 thermal expansion  alpha_xx  alpha_yy  alpha_zz 
    
 800.00  1.265114E05  1.240363E05  1.100857E05 
 1100.00  1.430439E05  1.403192E05  1.200152E05 
    
 mean value  Dal_xx  Dal_yy  Dal_zz 
    
 1.202111E05  5.24  3.18  8.42 
 1.344594E05  6.38  4.36  10.74 

    
 shearextension coeff eta_i,12  eta_i,13  eta_i,23 
    
 800.00 
 eta_1j  5.340722E08  5.363152E09  4.097442E02 
 eta_2j  3.741822E08  4.137568E09  4.043315E02 
 eta_3j  8.072154E09  4.344513E10  2.375421E04 
 1100.00 
 eta_1j  6.356865E08  9.550613E09  4.018362E02 
 eta_2j  4.509562E08  8.369636E10  3.952047E02 
 eta_3j  7.773691E09  4.573540E08  2.822128E04 
As expected, the equivalent composite material presents, in a first approximation, a transverse isotropic behaviour with a pronounced stiffness increase in fibre direction. Otherwise, the effective mechanical properties decrease with the temperature. As a thermoelastic homogenization analysis has been specified, HOMAT provides the effective thermal expansion coefficients in the structural directions. Due to the full anisotropy of the Hooke matrix, effective shearextension coefficients ^{14}, named \eta_{ij} with i \neq j, are given in the .hres
result file.
 In order to underline that the effective mechanical properties are expressed in the structural RVE axis system (x, y, z), each engineering property (e.g. Young, shear, Poisson, thermal expansion and eigenstrain) is doted by subscripts of type
xx
,xy
,yy
, etc.  If an isotropic eigenstrain is specified per phase, like in solid phase transformation simulations with MICRESS, the generalized effective volumetric eigenstrain tensor \kappa_{ij} is provided in addition to the thermal expansion coefficients. In the present analysis
VOLAVRG
is set toNO
, otherwise the volume averaged expressions of all the thermal and mechanical properties would have been provided. The correspondingPROBLEM_NAME.mat
file for a further thermal and/or thermoelastic analysis at the macroscale with Abaqus is given in Ex. 27. The symmetrised Hooke matrix is provided via the material keyword*ELASTIC
. Please note also the large analogy between the present Abaqus definition and the HOMAT material definition in Sec. Definition of the Materials Properties. As the special output optionENG=YES
is activated, HOMAT writes also the effective orthotropic Young and shear moduli and Poisson coefficients for the homogenized material. The anisotropic Hooke matrix suffices ABAQUS. Thus, the effective engineering constants are not read by ABAQUS but might be useful for later postprocessing.
,\*,\*
,\*,\*
,\*,\* NUMBER OF SPECIFIED MATERIALS: 1
,\*,\*
,\*MATERIAL, name = RVE
,\*,\*
,\*DENSITY
4.742357E03, 800.00
4.727331E03, 1100.00
,\*CONDUCTIVITY, TYPE=ORTHO
3.434308E02, 3.434307E02, 3.937139E02, 800.00
3.133230E02, 3.133228E02, 3.769298E02, 1100.00
,\*,\*
,\*SPECIFIC HEAT
6.559209E01, 800.00
7.480665E01, 1100.00
,\*ELASTIC, TYPE=ANISOTROPIC
3.001139E+05, 1.328911E+05, 3.000802E+05, 1.212415E+05, 1.212062E+05, 3.190173E+05, 3.669287E03, 1.808659E03,
1.304807E04, 8.040503E+04, 4.801965E04, 2.751873E04, 5.436726E06, 1.892165E+03, 8.446847E+04, 2.628284E+03,
2.557895E+03, 7.569334E02, 3.239438E04, 4.669365E04, 8.550840E+04, 800.00
2.694108E+05, 1.192916E+05, 2.693842E+05, 1.084949E+05, 1.084875E+05, 2.938578E+05, 3.875977E03, 1.909514E03,
1.194593E04, 7.134687E+04, 5.284290E04, 1.202110E03, 4.097078E03, 1.782636E+03, 7.873757E+04, 2.393042E+03,
2.312515E+03, 9.735707E02, 1.753401E04, 2.586395E04, 7.945887E+04, 1100.00
,\*,\*
,\*,\* 2. orthotropic variant neglecting shearextension coupling
,\*,\*
,\*ELASTIC, TYPE=ENGINEERING CONSTANTS
2.242213E+05, 2.242187E+05, 2.511521E+05, 3.422293E01, 2.800512E01, 2.798895E01, 8.036740E+04, 8.442893E+04,
8.542898E+04, 800.00
2.018505E+05, 2.018355E+05, 2.332937E+05, 3.459106E01, 2.791244E01, 2.791171E01, 7.131272E+04, 7.869988E+04,
7.938645E+04, 1100.00
,\*,\*
,\*YOUNG MODULUS, TYPE=ORTHO
2.242213E+05, 2.242187E+05, 2.511521E+05, 800.00
2.018505E+05, 2.018355E+05, 2.332937E+05, 1100.00
,\*,\*
,\*SHEAR MODULUS, TYPE=ORTHO
8.036740E+04, 8.442893E+04, 8.542898E+04, 800.00
7.131272E+04, 7.869988E+04, 7.938645E+04, 1100.00
,\*,\*
,\*POISSON RATIO, TYPE=ORTHO
3.422293E01, 2.800512E01, 2.798895E01, 800.00
3.459106E01, 2.791244E01, 2.791171E01, 1100.00
,\*,\*
,\*EXPANSION, TYPE=ORTHO
1.265114E05, 1.240363E05, 1.100857E05, 800.00
1.430439E05, 1.403192E05, 1.200152E05, 1100.00
,\*,\*
2D homogenization results¶
As the effective material can be either isotropic, orthotropic, cubic or generally anisotropic and the homogenization analysis can be of 3D or 2D nature, there exists several different output formats of the ASCII .hres
file. In the present manual, only two 2D variants: results of plane strain or plane stress homogenization, are briefly outlined.

Plane strain homogenization
The plane strain assumption (\epsilon_{33} =0) leads to following strain and stress tensors:
\begin{equation} \label{eq:str_planestrain} \epsilon^{\text{t}} = \left( \epsilon_{\text{xx}},\epsilon_{\text{yy}},\ 0,\ \gamma_{\text{xy}} \right)\,,\ \ \ \ \ \ \sigma^{\text{t}} = (\sigma_{\text{xx}},\ \sigma_{\text{yy}},\ \sigma_{\text{zz}},\ \tau_{\text{xy}})\,. \end{equation}The corresponding symmetric Hooke matrix, in 4\times 4 notation, has the following expression:
\begin{equation} \label{eq:H__planestrain} \mathbb{H} = \begin{bmatrix} H_{11} & H_{12} & H_{13} & H_{14} \\ \bullet & H_{22} & H_{23} & H_{24} \\ \bullet & \bullet & H_{33} & H_{34} \\ \bullet & \bullet & \bullet & H_{44} \\ \end{bmatrix}\,. \end{equation}The output format in the
.hres
file takes the reduced size (4\times 4) of the Hooke matrix into account, as illustrated by Ex. 28 of the UD reinforced glassepoxy composite, whose properties are given in Sec. Definition of the Materials Properties. As the outplane shears \tau_{13} and \tau_{23} are both zero, only the shear modulus in the crosssection plane 12 (or xy) is printed in the.hres
file.    HOMOGENISATION RESULTS  EFFECTIVE PROPERTIES     SIMULATION: per3D_scal05_evz     number of layers : 1   number of homogenisations : 1           Hooke matrix 4*4  H_i1  H_i2  H_i3        H_1j  1.443622E+04  4.898715E+03  5.274505E+03   H_2j  4.898715E+03  1.443606E+04  5.274471E+03   H_3j  5.274505E+03  5.274471E+03  4.626249E+04   H_4j  1.344940E08  4.578791E07  5.370989E08        part II  H_i4          H_14  1.445253E08     H_24  1.148775E08     H_34  5.451318E09     H_44  3.065298E+03               Shear moduli  G_xy           3.065298E+03           Young moduli  E_x  E_y  E_z         1.250001E+04  1.249987E+04  4.338476E+04        mean value  DE_x  DE_y  DE_z        2.279488E+04  45.16  45.16  90.33         Poisson coefficients  nu_xy  nu_xz  nu_yz         3.106216E01  2.727963E01  2.727974E01        mean value  Dnu_xy  Dnu_xz  Dnu_yz      

Plane stress homogenization
In a similar way, the plane stress assumption leads to the following stress and strain tensors:
\begin{equation} \label{eq:str_planestress} \epsilon^{\text{t}} = \left( \epsilon_{\text{xx}},\epsilon_{\text{yy}},\ \ \gamma_{\text{xy}} \right)\,;\qquad\sigma^{\text{t}} = (\sigma_{\text{xx}},\ \sigma_{\text{yy}},\ \ \tau_{\text{xy}}) \end{equation}The corresponding anisotropic Hooke matrix, expressed in Voigt notation, has following expression:
\begin{equation} \label{eq:H_planestress} \mathbb{H} = \begin{bmatrix} H_{11} & H_{12} & H_{13} \\ \bullet & H_{22} & H_{23} \\ \bullet & \bullet & H_{33} \\ \end{bmatrix}\,. \end{equation}The output in the
.hres
file takes the reduced size of the Hooke matrix (3\times 3) into account, as illustrated by Ex. 29 of the UD reinforced glassepoxy composite. Hooke matrix 3*3  H_i1  H_i2  H_i3        H_1j  1.027775E+04  6.249302E+02  0.000000E+00   H_2j  6.248244E+02  1.027799E+04  0.000000E+00   H_3j  0.000000E+00  0.000000E+00  3.932720E+03         Young/shear moduli  E_x  E_y  G_xy         1.023976E+04  1.024000E+04  3.932720E+03        mean value  DE_x  DE_y         1.023988E+04  0.00  0.00          Elastic moduli  C_11  C_12  C_44         1.027787E+04  6.248773E+02  3.932720E+03        Eng. moduli  E  nu  Zener coeff.        1.023988E+04  6.079763E02  8.148189E01         Poisson coefficients  nu_xy  nu_yx          6.079763E02  6.079903E02  
As in this analysis the parameter
CUBVAL
of the keyword/PARA
is set to0.1
, the transverse quasiisotropy of the UD reinforced composite leads HOMAT to provide also the effective elastic properties of the equivalent cubic material. Indeed, the three elastic constants C_{11}, C_{12} and C_{44} are given as well as the mean Young modulus, the mean Poisson coefficient and the Zener anisotropy coefficient, defined by\begin{equation*} Z = \frac{2C_{44}}{C_{11}  \ C_{12}}\,. \end{equation*}Moreover, in this run
VOLAVRG=YES
, therefore the volume averaged Hooke matrix is also printed. For the thermal expansion coefficient, not only the volume averaged expression is evaluated via the thermal stresses and printed but also the value obtained by the direct application of the mixture rule to the thermal expansion coefficient without using the thermal stress expression. The corresponding results are shown in Example 30 below. th. exp. vol. aver.  alpha_11  alpha_22          4.406312E06  4.406312E06          th. expan. mixture  alpha_11  alpha_22          9.587308E06  9.587308E06  
Note: The different output of both 2D homogenization analyses with some different parameters (
CUBVAL
,VOLAVRG
) shows the great variety and the complexity of the output file.hres
.
Results in a local axis system¶
Since version 6.0 the user has the possibility to specify a local axis system in which the calculated effective Hooke matrix is printed in the PROBLEM_NAME.hres
and PROBLEM_NAME.mat
files. To illustrate this new output a more industrial example of a hollow wireframe sandwich, built of two Nicrofer plates which surround woven firing and chaining wires (see Fig. 17), is adopted here. This hollow wire frame is placed on the stator of a steam turbine and it is cooled by steam at 425°C. The thermal and mechanical homogenizations are performed at 650°C. Due to the ECAP drawing process, the wires present elongated grains in the extrusion direction. Therefore, their mechanical properties are orthotropic. For each wire a local axis system in which the material properties are specified via the coordinates of two nodes (see Ex. 31), is used. In this example, the local axis system of the chain wire is chosen via parameter OUT_ORIE
on the /ORIENTATION
keyword for additional output of the effective anisotropic Hooke matrix of the cooled hollow multilayer. Note that this additional output requires the setting of LOC_AXIS=YES
as parameter of the /OUTPUT
keyword.
#
# Output definition
#
/OUTPUT, FORMAT=VTK, TYPE=ASCII, IOP=2, VOLAVRG=YES, VISLAV=YES,
LOC_AXIS=YES, SHEAR_TENS=YES
U
/ORIENTATION, NAME=wire_X_axis, DEFINITION=COORDINATES, OUT_ORIE
0., 2.0, 0., 0., 0., 2.0
/ORIENTATION, NAME=wire_Y_axis, DEFINITION=COORDINATES
1.0, 0.0, 0.0, 0.0, 0.0, 1.0
#
# Definition of the solid sections
# per solid section, specify the desired analysis type
#
/SOLID SECTION, NUMBER=1, ELSET=SCHUSSDRAHT, NAME=firing_wire,
MATERIAL=HT6025,
ORIENTATION=wire_Y_axis, TEMP=NO, MECH=HOM
#
/SOLID SECTION, NUMBER=2, ELSET=KETTDRAHT, NAME=chain_wire,
MATERIAL=HT6025,
ORIENTATION=wire_X_axis, TEMP=NO, MECH=HOM
In fact, the local axis 1 is directed along structural axis Y; whereas local axis 2 corresponds to the structural direction Z and the 3^{d} orthogonal axis is along structural axis X. Therefore, the local Young, shear moduli, listed in Ex. 32 correspond to a simple permutation of the structural values. This permutation can easily be verified by the user.

Output in the structural axis system:
 Output in the local axis system Shear moduli  G_xy  G_xz  G_yz         4.255631E+04  6.807693E+04  6.807691E+04        mean value  DG_xy  DG_xz  DG_yz        5.957005E+04  28.56  14.28  14.28         Young moduli  E_x  E_y  E_z         1.174384E+05  1.176799E+05  1.607221E+05        mean value  DE_x  DE_y  DE_z        1.319468E+05  11.00  10.81  21.81         Poisson coefficients  nu_xy  nu_xz  nu_yz         2.787595E01  4.118723E01  3.931044E01        mean value  Dnu_xy  Dnu_xz  Dnu_yz        3.612454E01  22.83  14.01  8.82 
 local shear moduli  G_12  G_13  G_23         6.807691E+04  4.255631E+04  6.807693E+04        mean value  DG_12  DG_13  DG_23        5.957005E+04  14.28  28.56  14.28         local Young moduli  E_1  E_2  E_3         1.176068E+05  1.605981E+05  1.174333E+05        mean value  DE_1  DE_2  DE_3        1.318794E+05  10.82  21.78  10.95         local Poisson coeff.  nu_12  nu_13  nu_23         2.879694E01  2.789950E01  4.116669E01        mean value  Dnu_12  Dnu_13  Dnu_23        3.262104E01  11.72  14.47  26.20 
Stokes flow homogenization input and results¶
For each specified layer or classically for the whole multilayer, HOMAT provides the calculated Darcy permeability tensor in the structural axis system of the RVE and usually in mm^{2}. Additionally, the permeability tensor is also written in the local axis system, if such a local axis system has been specified via parameter LOC_AXIS
on the /OUTPUT
keyword and via parameter OUT_ORIE
to select its orientation. As the definition of a Stokes flow homogenization problem is new since version 6.0, not only the permeability results are outlined for a benchmark test, but also the specific commands of the problem definition. The flow around and parallel to a cylindrical fiber of radius 2 mm in a cubic box of dimension 10 \times 10 \times 10 mm^{3} is taken here as example. Ex. 34 outlines the specific commands for Stokes flow homogenization. Note that at least one domain has to be defined as fluid domain via following allowable solid section names: FLUID
, AIR
, GAS
, MELT
, LIQUID
, POROUS
or by following German names LUFT
and SCHMELZE
. Note: No other names are accepted.
#
/PROBLEM, TYPE=STOKES, HOM_MOD=0, PER_MOD=3
#
/OUTPUT, FORMAT=VTK, TYPE=ASCII, IOP=2, VOLAVRG=YES, ENG=NO, VISLAV=YES, LOC_AXIS=NO
V, P
#
/INCLUDE, NAME=FibBL_r2_bubLi.hgeo.gz
#
MATERIAL, NAME=fiber, TYPE=ISO
/INCLUDE, NAME=silica.hmat
#
/MATERIAL, NAME=AIR, TYPE=ISO
/INCLUDE, NAME=luft_mm.hmat
#
/SOLID SECTION, NUMBER=1, ELSET=SetMatrix, NAME=fluid, MATERIAL=AIR,
TEMP=NO, MECH=NO, FLOW=CALC
#
/SOLID SECTION, NUMBER=2, ELSET=SetFiber, NAME=FIBER, MATERIAL=fiber
TEMP=NO, MECH=NO, FLOW=NO
#
In order to respect the general definition of the /SOLID SECTION
keyword, also the material has to be specified, nevertheless the defined properties (density, viscosity, thermal conductivities, ...) are not used by HOMAT. HOMAT solves three (in 3D) or two (in 2D) special Stokes flow problems, where the density and the viscosity are both taken as unitary and a unitary volumetric boost force is applied in the structural directions X, Y and Z, respectively. Independent of the user specification FLOW=YES
or NO
for a solid domain, HOMAT exclude the solid domains from the reduced system, except the interface nodes, where noslip boundary conditions are always applied. In Example 35 the predicted Darcy permeability tensor for the cylindrical fiber RVE is outlined. HOMAT 6.0 provides the principal permeabilities and the corresponding eigenvectors. Moreover, it computes the R anisotropy coefficient defined by the principal permeabilities K_{1} > K_{2} > K_{3} via the expression ^{19}

 
 HOMOGENISATION RESULTS  EFFECTIVE PROPERTIES 
 
 SIMULATION: FibBL_r2_bubLi 
 porosity: 87.45 % 
 
 number of layers : 1 
 number of homogenisations : 1 
 

    
 permeability tensor  K_X  K_Y  K_Z 
    
 comp. X  5.399124E+00  4.867793E12  1.735329E15 
 comp. Y  4.123367E13  5.399125E+00  6.511195E16 
 comp. Z  3.288843E10  2.218531E10  6.668920E+00 
    

    
 princ. perm. R_ani  K_1  K_2  K_3 
    
 8.99775E01  6.668920E+00  5.399125E+00  5.399124E+00 

    
 unit vector  K_1  E1_1  E2_1  E3_1 
    
  1.295035E10  8.735858E11  1.000000E+00 
    
 unit vector  K_2  E1_2  E2_2  E3_2 
    
  2.306980E06  1.000000E+00  8.735806E11 
    
 unit vector  K_3  E1_3  E2_3  E3_3 
    
  1.000000E+00  2.306980E06  1.295037E10 
Note that HOMAT also provides the porosity of the RVE as additional information. The predicted effective Darcy permeabilities are, as expected, transverse isotropic and the offdiagonal terms are nearly zero. The permeability parallel to the fibers, denoted K_{\parallel}, is larger than the transverse one due to less flow resistance parallel to the fibers, characterized by a lower tortuosity. These predictions are achieved with a coarse hexahedral mesh of 37120 elements and 41517 nodes but with a boundary layer around the fiber (see Figs. 20 and 21 of the pressure distribution around the central fiber). The permeability K_{\parallel} is in excellent agreement with the analytical solution of Tamayol and Bahrami ^{20} for the flow parallel to a square arrangement of cylindrical fibers: K_{\parallel}^{a} = 6.6731 mm^{2;}; whereas the transverse permeability prediction P_{\perp}= 5.3991 mm^{2} overestimates the model value K_{\perp}^{a}= 3.696 mm^{2} of Tamayol and Bahrami ^{21}, leading to a less pronounced orthotropic permeability behaviour, in agreement with experimental results at high porosities (here \epsilon = 0.874).
Note: If the parameter LOC_AXIS=YES
is specified, the permeability tensor is also printed in the local axis system (1, 2, 3). Compared to the output of Ex. 35, the descriptor becomes now local perm. tensor instead permeability tensor and the permeability components are noted now K_1
, K_2
and K_3
instead K_X
, K_Y
and K_Z
, indicating thus the output in the local axis system (1, 2, 3).
Field Results¶
In this section the previous example of the sapphire fibre in a metal matrix is reconsidered to detail the content of the different field output files in the case of a combined thermal and thermoelastic homogenization. Field results are characterized by nodal values, defined for each considered layer. They describe the local variation of the field on the layer or RVE and are written in the VTK format. They can be visualized via the opensource program Paraview.
Thermal field results¶
If the parameter T
or ALL
is specified in the /OUTPUT
keyword, each thermal mono or multilayer homogenization produces two or three different VTK files.
 The 1^{st} file contains all the generated implicit fluxes at the domain interfaces, where a jump of material properties exists. If the thermal conductivity varies inside each domain with the temperature, it contains also their volumetric contribution. The number of flux terms depends on the homogenization type: either 3 in the 3D case or 2 in the 2D case. The corresponding keyword in the VTK file is
IMPL_FLUX
. Since version 5.2, the implicit fluxes are saved as vector solutions in the VTK file. Its nomenclature is the following: to the problem name,fas30_fine
, is added extension_flux
. Then, two numbers are added: the 1^{st} one specifies the considered layer and the second one, the considered homogenization run. Each set temperature on the RVE (keyword/TEMP DOMAIN
) leads to a separate homogenization analysis. In the present case, as two temperatures are specified (T=800°C and T=1100°C) on one layer, HOMAT writes following two implicit flux files for visualization purpose:fas30_fine_flux11.vtk
andfas30_fine_flux12.vtk
: the 1^{st} corresponds to the thermal homogenization at 800°C and the 2^{d} one to the homogenization at 1100°C. 
The second file contains the periodic microscopic displacements x^{r} , with r=1, 3 in the 3D case or r=1, 2 in the 2D case on the undeformed RVE. Each displacement component is the solution of a linearized thermal system. The corresponding keyword in the VTK file for this vector solution is:
DISPL_TH
. Its nomenclature is similar to the one for implicit fluxes: HOMAT adds to the problem name the extension_hola
followed by two numbers: the layer and homogenization indicators. Thus, for the metal matrix composite example, HOMAT generates the following two files with microscopic displacements:fas30_fine_hola11.vtk
andfas30_fine_hola12.vtk.
Note that these vtk files have a common part composed by a header, the definition of the undeformed finite element mesh and the local nodal values of thematerial_ID
parameter. This parameter allows the drawing of the different domains or grains of the RVE. In Fig. 18 (left) the implicit flux in X direction is drawn for the UD reinforced metal matrix composite and in Fig. 18 (right) the corresponding microscopic displacement. Due to the higher thermal conductivity of the FG75 matrix than the sapphire fiber, axial implicit fluxes are generated in X direction, which tends the fiber to be compressed along this axis. As the implicit fluxes are usually small, they are multiplied by a constant factor of 10^{3} in the flux result files.
Mechanical field results¶
If the parameter U
or ALL
is specified in the /OUTPUT
keyword, each thermoelastic mono or multilayer homogenization analysis produces the output of an implicit forces file and of an microscopic displacement field file on the undeformed RVE.
 The 1^{st} file contains all the generated implicit forces at the domain interfaces, characterized by a jump in the material properties. If the temperature varies inside each domain, it contains also a volumetric contribution (see Appendix Homogenization of Materials with Periodic Microstructure). Its nomenclature is similar to the implicit flux file, except, of course, its extension differs:
_hfor
instead of_flux
. In the present example, as two temperatures are specified (T=800°C and T=1100°C) on one layer, HOMAT writes following two implicit force files:fas30_fine_hfor11.vtk
andfas30_fine_hfor12.vtk
. The 1^{st} corresponds to the thermoelastic homogenization at 800°C and the 2^{d} one to the homogenization at 1100°C.  The second file contains again the periodic microscopic displacement fields \zeta^{rs} with \text{rs}=1, 6 in the 3D case, on the undeformed RVE. This file is characterized by the extension
_hoel
. Otherwise, its nomenclature is identical to the implicit forces file. Therefore, following two periodic displacement files are written:fas30_fine_hoel11.vtk
andfas30_fine_hoel12.vtk
.
The number and nomenclature of different vector solutions in each of these files depends not only of the nature of the homogenization analysis (pure mechanical or thermoelastic) but also on the homogenization model type: 3D or 2D (plane strain or plane stress model):
 3D: in the case of a pure mechanical homogenization, we have six forces and microscopic displacements, corresponding each to a unitary macrostrain: either E_{XX}, E_{YY}, E_{ZZ}, E_{XY}, E_{XZ} or E_{YZ}. Each force/displacement vector solution has 3 components in the VTK file. The implicit forces have the identifier
IMPL_FORC_
followed by the extension of the nominal strain (e.g.XX
). For example, HOMAT defines in the force filefas30_fine_hfor11.vtk
, the solutionsIMPL_FORC_XY
andIMPL_FORC_YY
, which means that an implicit force vector of 3 components is generated by the unitary strains E_{XY} and E_{YY} respectively. For the microscopic displacements on the undeformed RVE, the nomenclature is similar except, of course, the identifier name becomesDISPL_
. Therefore, in the filefas30_fine_hoel11.vtk
, following solution identifiers of microscopic displacements are then specified:DISPL_XY
,DISPL_YY
. In the case of thermoelastic homogenization, an additional implicit force solution due to the thermal expansion coefficient and the corresponding microscopic displacement solution are saved in the respective VTK files. They are characterized by specific identifier names:TH_EXP_FORCE
for the implicit force andDISP_TH
for the microscopic displacement field. Moreover, if volumetric eigenstrains are active, a second additional implicit force solution, with identifierTH_EIG_FORCE
is written in the VTK filePROBLEM_NAME_hfor.vtk
and the corresponding microscopic displacement field, with identifierDISP_EIG
is written in the VTK filePROBLEM_NAME_hoel.vtk
. In this special thermoelastic case, eight vector solutions are output.  2D plane strain model: in the case of a purely mechanical homogenization, only four implicit forces, namely F_{XX}, F_{YY}, F_{ZZ} and F_{XY} are evaluated and the two planar components X and Y are saved in the implicit force file with the same identifier
IMPL_FORC
as in 3D. In presence of a thermoelastic homogenization, one additional implicit force vector with two components, namedTH_EXP_FORC
and corresponding to the thermal expansion effect, is added. Moreover, if volumetric eigenstrains are active, a second additional implicit force vector solution, namedTH_EIG_FORCE
is introduced. The corresponding microscopic displacement file contents then either four (pure mechanical case), five (thermoelastic case, only thermal expansion coefficient) or six solutions (thermoelastic case with eigenstrain and th. expan. coeff.). The identifiers of these solutions are the same as for the 3D case, namely:DISPL_XX
,DISP_YY
,DISP_ZZ
,DISP_XY
,DISPL_TH
andDISPL_EIG
. 
2D plane stress model: three implicit forces, namely F_{XX}, F_{YY} and F_{XY} are expressed in a pure mechanical homogenization analysis. Again only the planar components are saved per force or displacement, leading to three stored solutions with 2 components per solution in the VTK file. The nomenclature and the identifier names are the same as for the 3D case. In a thermoelastic homogenization analysis, again additional inplane solutions are expressed for the thermal expansion coefficient and, if specified, for the volumetric eigenstrain in each VTK file, leading to a total number of either four or five vector solutions. In Fig. 19 (left) the microscopic displacement in X direction, induced by the unitary strain E_{XX} is drawn on the deformed RVE of the UD reinforced metal matrix composite and in Fig. 19 (right) the microscopic displacement in Y direction, induced by the unitary strain E_{YY} is outlined. As the elasticity of the sapphire fibre is stronger than the elasticity of the FG75 NiAl matrix, the implicit forces in X direction, induced by the macrostrain E_{XX}, tends to extend the fibre in X direction (formation of an ellipse). Otherwise, the implicit force in Y direction under macrostrain E_{YY} tends also to distort the fibre elliptically with the larger semiaxis in Y direction.
Remarks:

Due to the smallness of the thermal expansion coefficients, the corresponding microscopic displacements are also small and in order to better visualize these displacements they are scaled by a factor 10^{3}. Note that the volumetric eigenstrains are usually larger, so that their corresponding displacements are not scaled by default.

In order to represent a 2D vector solution, the user must create an additional solution by using the calculator option in Paraview. For example, for the microscopic displacement solution
DISPL_XX
, the user introduces in the calculator following expression:iHat*DISPL_XX_X + jHat*DISPL_XX_Y + kHat*0
This additional solution permits to use the “wrap by vector” option in order to visualize the desired planar vector solution, for example,
DISP_XX
. 
Since version 6.0, the user can also specify the parameters
P
andE
in the/OUTPUT
keyword. The parameterP
, corresponding to the local pressure field, is set, by default, if stabilized mixed finite elements are defined in the thermoelastic homogenization 3D or 2D analysis to model quasiincompressible material behaviour. The additional output parameterE
allows to display the periodic microstrains \vec{\varepsilon}(\vec{\varsigma}^{\text{rs}}) induced by the variation of the microscopic displacements \vec{\varsigma}^{\text{rs}} over the RVE. This way, the user can more easily analyse the evaluation of the effective Hooke matrix via expression \eqref{eq:eff_hoke_mat} and the effective eigenstress via expression \eqref{eq:eff_eigenstress}. 
In the presence of mixed finite elements, a scalar pressure field is saved for each resolved linear problem: The number of saved pressure fields in quasiincompressible materials varies from 6 to 8 in 3D homogenization, from 4 to 6 in 2D plane strain or from 3 to 5 in 2D plane stress homogenization. The file of periodic pressure fields is characterized by the extension:
_hdru
. Otherwise, its definition is identical to the nomenclature of the microscopic displacement file. For the example of the crystalline/amorphous lamella design at the nanoscale (see Fig. 9 in ^{2}), the following pressure file is written:bubN_draw80_hdru11.vtk
. 
As a purely mechanical homogenization is performed in 3D in this example, this file contains the following scalar pressure solutions, named:
PRES_XX
,PRES_YY
,PRES_ZZ
,PRES_XY
,PRES_XZ
, andPRES_YZ
. Moreover, as the parameterALL
was specified in the command file also the microscopic strain tensors are written for visualization purposes. The file of microscopic strains is characterized by the extension_hdefo
and in the present example it isbubN_draw80_hdefo11.vtk
. This file contains in 3D the six microstrain tensor compnentsSTRAIN_XX
,STRAIN_YY
,STRAIN_ZZ
,STRAIN_XY
,STRAIN_XZ
andSTRAIN_YZ
. In the case of thermoelastic homogenization, one or two additional strain tensors are written in thePROBLEM_NAME_HDEFO.vtk
file. They are characterized by the identifiersSTRAIN_TH
andSTRAIN_EIG
respectively.
Note: By analogy to ABAQUS, the microscopic strain tensors contain engineering shear strains with \gamma_{xy} =2 \varepsilon_{xy}. Paraview easily interprets the six components of each microscopic strain tensor as tensorial components. However, for plane strain homogenization analysis, Paraview denotes the four strain components as 0, 1, 2, 3 which leads to following correspondence: 0 = XX, 1 = YY, 2 = ZZ, 3 = XY. If 2D plane stress analysis is adopted, Paraview denotes the three stress/strain components X, Y and Z. In this case, we have simply following correspondence: X = XX, Y = YY and Z = XY.
Stokes flow field results¶
By default, for a fluid flow homogenization analysis the parameters V
and P
are specified on the /OUTPUT
keyword and HOMAT v6.0 prints for each mono or multilayer homogenization three different VTK files:
 The first file contains all the generated implicit volumetric boost forces, corresponding to an unitary pressure difference applied in the structural direction X, Y, and Z respectively. The number of boost terms depends on the homogenization type: either 3 in the 3D case or 2 in the 2D case. The corresponding keyword in the VTK file is
BOOST_B_*
. The star*
stands either forX
,Y
orZ
in 3D, specifying the structural direction of the boost forces. These boost forces are save as vector solutions in the corresponding VTK file. Its nomenclature is the following: to the problem name,FibBL_r2_bubLi
, is added extension_boost11
. To preserve the same definition as for the VTK definition of thermal and mechanical homogenizations, two numbers are also added: the first one specifies the considered layer and the second one, the considered homogenization. Note that they are always 1 in the case of Stokes flow homogenization.  The second file contains the periodic microscopic velocities \omega^{r}, with r=1, 3 in the 3D case or r=1,2 in the 2D in fluid domains of the RVE. Of course, these velocities are zero in the solid domains of the RVE. Each velocity vector has three (or two) components. They are the solution of the special Stokes problems. HOMAT adds to the problem name the extension
_hvel11
to specify the corresponding VTK file. The three or two velocity vector solutions in this VTK file have as identifier:VELO_B_*
, where again*
stands forX
,Y
orZ
, respectively. Note thatB
in the identifier recalls that boost forces in the respective direction X, Y or Z are responsible for these velocities. 
The third VTK file contains the periodic pressure fields in the fluid domains of the heterogeneous material. A scalar pressure field is saved for each resolved Stokes flow problem. The number of saved pressure fields is either 3 in 3D or 2 in 2D fluid flow homogenization. The file of periodic pressure fields is characterized by the extension:
_hdru11.
In this VTK file, the different scalar pressure fields are denoted via identifierPRES_XX
,PRESS_YY
andPRESS_ZZ
in 3D. For the example of the flow transverse to a cylindrical fiber of radius 2 mm in a cubic box (10*10*10 mm^{3}) the pressure fieldsPRES_XX
andPRES_YY
are drawn in Figs. 20 and 21, respectively.
Localization Results¶
Microstrains and stresses are evaluated on the RVE in selected, critical points of the thermomechanical structural analysis of the sample at the macroscale. For each of these locations, total microstrains are evaluated via expression \eqref{eq:loc_str} in the Gauss points of the elements and elastic microstresses are calculated via expression \eqref{eq:loc_str_1} in the same Gauss points. Independent of the finite element discretization, a linear extrapolation scheme and a dedicated averaging technique are used to express their corner node values. As strains are continuous fields through material interfaces, all the nodal contributions on both interface sides are taken into account in the weight function; whereas for discontinuous fields, like stresses, only contributions from the same material domain are considered.
For each localization run, the total strain tensor, an equivalent strain as well as the elastic stress tensor and the equivalent von Mises stress are saved in the element nodes of the RVE or of the selected layer. By specifying the monolayer option using the keyword /DEF LAYER
, it is possible to select the layers for which microstresses and strains are evaluated. This feature increases the efficiency of the program, as less CPU time is required to gain the desired microstructure results.
Note that, by analogy to ABAQUS, the total strain tensor contains engineering shear strains with \gamma_{xy} =2 \epsilon_{xy}.
The number of stress/strain components depends on the localization analysis type: six in 3D plus the equivalent values. In 2D plane strain analysis they are four plus the equivalent term and in 2D plane stress analysis, three plus the equivalent stress/strain. In the VTK files appear either the identifier STRAIN
or STRESS
for both tensors and STRAIN_EQ
or STRESS_EQ
for both scalar equivalent measures. Note that in 3D, Paraview interprets easily the six components of both tensors as tensorial components. But, for a 2D plane strain homogenization analysis, Paraview denotes the four stress/strain components as 0, 1, 2, 3, which leads to following correspondence:
The microstrain and stress files have the following nomenclature: after the problem name, the extension _loeps
specifies a microstrain VTK file, which contains the components of the strain tensor as well as the equivalent strain; the extension _lostr
specifies a microstress VTK file, which contains the component of the stress tensor, the equivalent stress and, for mixed elements, the pressure. To this extension is then added the number of the localization analysis. In the header of each VTK file, HOMAT writes the selected critical point of the macrosimulation (Gauss point of the considered macroelement, part of specified element set).
As an example, consider the localization analysis performed after the UO forming simulation of a ferritepearlite pipeline tube (see Fig. 1). In the critical area ^{22} two Gauss points nearest the upper and lower surface of the X65 steel sheet are selected, where total microstrains and stresses are expressed in the ferritepearlite microstructure via the localization procedure. The VTK file names are respectively: rve2_loeps_1.vtk
, rve2_loeps_2.vtk
for the microstrains and rve2_lostr_1.vtk
, rve2_lostr_2.vtk
for the microstresses.
In Fig. 22 (left) the most significant shear strain component is drawn on the RVE and in Fig. 22 (right) the equivalent stress distribution in the pearlite inclusions is illustrated. This strain distribution shows clearly that some elongated pearlite inclusions, oriented in rolling direction, are strongly sheared. As pure elastic localization analysis is performed and no plasticity is taken into account, the equivalent stresses in the pearlite inclusions are overestimated.
Log File Output¶
In this paragraph a brief overview of the most relevant output on the screen or  if redirected  in the log file, named either PROBLEM_NAME.hout
, for a homogenization analysis, or PROBLEM_NAME.lout
for a localization run, will be presented briefly. This output is so variable and depends on all possible combinations of input data, that a rigorous, detailed description of the log file output is impossible. Here some characteristic messages and information are highlighted. These messages will help the user to verify the correctness of its analysis. First, some common messages are discussed, then specific messages of a thermal homogenization run are presented, followed by messages of a thermoelastic homogenization analysis and of the Stokes flow homogenization and, finally, by the description of the specific output of a localization run.
Common output¶
Each HOMAT run starts by the interpretation of the command file. The corresponding log file messages are selfexplanatory. The first important check is performed in HLIGEO
to verify if the bounding box around the RVE is correctly and that its dimensions in LENCHAR
are correct expressed in the adopted unit length. The determination of nodes belonging to a periodic surface pair is a delicate operation, mainly in presence of different finite element (FE) meshes on both linked periodic surfaces. The user can check at first the correctness of the mean and minimum element wedge length used in the projection algorithm of the slave nodes on the master surface. Then, the output of routine PROPER
for each periodic surface is important, as it gives the number of directly corresponding nodes and the number of nodes, which require, after projection on the master surface, an interpolation within the element containing the projection point. Globally, we can say that the more nodes need to be interpolated, the slower is the convergence of the iterative solver. Therefore, in order to reduce the number of nodes to be interpolated for incompatible surface meshes, the user can specify a larger value of the EPSCOORD
parameter (see keyword /NUMERICS
), for example 5.E02
instead of 1.E02
. Sometimes, it can be interesting to check the number of detected double nodes per wedge of the RVE, provided by routine HLIECK
. For the dimension of the homogenization or localization problem, the total number of Gauss points provided by routine FLAEPER
is an important quantity. If only one layer is specified in the model (either multi or monolayer homogenization), the volume and its fraction are expressed for each phase or material. These results provide important information as they measure, for example, the discrepancy between the discretized phases in the RVE/layer and the experimental or theoretical volume fraction of each phase, mainly if the MICRESS file .korn
file is used to separate the individual grains of the microstructure. Independent of the considered homogenization analysis, the total surface area of each periodic surface is printed and help the user to verify the model's correctness. For each nonzero righthand side vector, the reduced global linear system is solved iteratively. The output frequency of intermediate solver results, like the value of the convergence criterion and the actual norm of residuum, is governed by the OUTFREQ
parameter of keyword /NUMERICS
. The user can check via these outputs the convergence progress. The value of residuum and relative convergence criterion at the best iteration has to be analysed carefully. A too large relative convergence criterion (\epsilon > 0.01 ) is an indication of a badly configured global system. Bad systems occurs mainly in presence of large implicit forces or fluxes due to strong material discontinuity (e.g. crystalline/amorphous lamella) and sometimes in presence of 2D periodicity inducing large microscopic fluctuation displacements on the RVE: Please, check if a 3D periodicity B.C. cannot help to circumvent these too large fluctuation displacements on the RVE.
Thermal homogenization output¶
The norm of the implicit source terms is expressed at the end of the element generation in the routine GENVOTH
. If a constant temperature pro domain is specified, no source terms are generated and their values have to be zero. Then, only implicit forces at the boundary between two different materials or phase are generated. These fluxes have to be in equilibrium in each structural direction. Please control their equilibrium. For each domain, the volume and the mean temperature are provided. The size of the reduced system, defined by the number of Degrees Of Freedom (DOF) is given by the routine RANDIR
, please check this relevant information for the memory allocation. For each structural direction the maximum and minimum value of the implicit fluxes and their nodal location are outlined. If the normal at the material interfaces has no component in a specific direction, for example, in the third one (case of UD reinforced fiber composite), then the corresponding min. and max. values are zero and no contribution to the total implicit forces is expressed. Similarly, after the resolutions of the global reduced system, the min. and max. values of the predicted microscopic displacement field are written for each component. Eventually, the effective thermal conductivity tensor and the other effective thermophysical properties (thermal capacity, specific heat, density) are provided. The corresponding volume averaging values are also printed.
Mechanical homogenization output¶
The mechanical finite element generation starts with the output of the Hooke matrix in the structural axis system of the 1^{st} element of each domain. Please carefully verify these Hooke matrices. Again implicit volume forces are only evaluated in presence of temperature variations within a RVE domain or when the spherulite model is used. In this case, different lamella orientations in each Gauss point induce implicit volume forces. Then, for each identified interface with a jump of mechanical properties, the evaluated implicit mechanical forces are given for each component. In the case of thermoelastic homogenization, also the total implicit forces due to unit thermal stresses, are printed. For each structural direction the maximum and minimum value of the implicit forces and their nodal location are given for the mechanical homogenization. Please check if the number of solved microscopic problems corresponds to the desired one for the defined homogenization analysis. Please check the number of DOF of the reduced system, printed by routine RANFIX
; it directly impacts the memory allocation for the considered homogenization run. After the resolution of the different global systems, the min. and max. values of each predicted microscopic displacement component are given. Subsequently, the effective Hooke matrix before symmetrisation, the symmetric Hooke matrix and its inverse, the symmetric compliance matrix, from which the effective engineering moduli are derived, are provided. Thereafter, the Hooke matrix in the specified local axis system and/or the volume averaging engineering moduli are printed.
Stokes flow homogenization output¶
For each fluid domain, the total volumetric boost forces for each structural direction are provided by the GENFLOW
routine. In addition, for each domain (solid or fluid) its volume is written. Further, for each material or phase, the volume fraction is printed, specifying the porosity of the porous medium. Please check this porosity for its correctness. After redistribution of the master nodal values, via parameter VISLAV=YES
, the max. and min. boost forces are provided. Due to nonslip boundary conditions at fluidsolid interfaces, the reduced system formed of velocity and pressure DOF in the fluid domains is further reduced and the corresponding total number of DOF of the Stokes flow problem is printed. Then, after resolution of the three or two special Stokes problems, the min. and max. values of each predicted microscopic velocity solution and of the microscopic pressure field are written. Finally, the predicted Darcy permeability tensor in the structural axis system is written in the log file.
Localization output¶
The begin of the localization output is identical with the log output of the mechanical homogenization run, as the microscopic displacement fields have to be determined for the calculation of the microscopic strains and stresses respectively. After the resolution of the different microscopic problems, starts the real localization analysis for the different specified macrostrain states. For each localization analysis, the mean internal work per layer or for the whole multilayer is evaluated and printed. In the case, if the internal elastic work is expressed for the complete RVE, its value can be compared with the deformation work at the macroscale and can thus be used to validate the HillMandel principle ^{23}. Moreover, the mean internal work for each phase or material is given, showing so the contribution of each phase/material to the internal work of the layer or the RVE. Finally, the routine LOCWRITE
performs a similar output of min. and max. nodal values, for each component, of the total microstrain field and microstress field on the RVE or on the considered layer.
Transfer of macroscale results from ABAQUS to HOMAT¶
In order to perform a linear elastic localization analysis with HOMAT, the user has to provide via the keyword /INPUT FILES
, the name of the ASCII file of the linked thermomechanical macrosimulation of the part with ABAQUS. This ASCII file has always the suffix .txt
. A dedicated python script, named Aba2Hom.py
, has been written in order to extract from the data base of macrosimulation (.odb
) field results (strains, stresses, temperature, von Mises stresses,...) in critical regions of the sample at selected load steps and iterations. Please look at the corresponding user guide of this python script to perform the desired extraction in critical regions of the specimen.
Here, the focus is given on the description of the transfer file NAME.txt
, which is read by HOMAT. Example 36 shows such a file extracted from a plane strain thermoelastic analysis of a 10° degree sector of a cooled cylindrical combustion chamber wall subjected to a stationary thermal gradient (see Fig. 23, left). The final deformation and the radial strain are plotted in Fig. 23 (right). The region of interest is magnified in Fig. 24. Strain, stress and temperature results are extracted by applying Aba2Hom.py
with the analysis type CritSelect
. for the elements, labeled 10, 30, 50, 79 and 80. The generated transfer file is shown in Example 36.
/HEADER
macrostrains extracted in the center of combustion chamber wall
Step1
ODB Name:orth_lin.odb
ABAQUS Version 2019
HOMAT Version 6.0
Thu April 2 10:23:30 2020
/PARA, ABA_MOD=1, NEL=5, NPG=4
/ELSET=CRIT_TBC+BOND
/STRAIN
10,1,0.00033,0.01487,0.00000,0.00013
10,2,0.00033,0.01487,0.00000,0.00003
10,3,0.00047,0.01486,0.00000,0.00014
10,4,0.00047,0.01486,0.00000,0.00004
......................................
79,1,0.01610,0.01481,0.00000,0.00001
79,2,0.01617,0.01481,0.00000,0.00001
79,3,0.01610,0.01481,0.00000,0.00000
79,4,0.01617,0.01481,0.00000,0.00000
.....................................
/TEMP
10,1,931.26550293
10,2,931.26550293
10,3,907.40246582
10,4,907.40246582
..................
Aba2Hom.py
six header lines. In the first two lines the user can specify a title describing the performed multiscale investigation. On the third line, the selected ABAQUS analysis step is given, whereas, in 4^{th} line the name of the extracted ABAQUS data base is written, here orth_lin.odb
. Eventually, in line 5^{th} and 6^{th} the used program version of ABAQUS and HOMAT are provided.
The keyword /PARA
specifies the used ABAQUS discretization model, via the mandatory model parameter ABA_MOD:
ABA_MOD=0
: 3D macroscale analysis with six strain componentsABA_MOD=1
: 2D plane strain, axisymmetric or plane stress analysis with four strain componentsABA_MOD=2
: 3D membrane or shell analysis with three strain componentsABA_MOD=3
: axisymmetric shell elements having two strain componentsABA_MOD=4
: truss elements having only one axial strain component
This keyword also provides the number of selected elements via the parameter NEL
and the number of Gauß points per element via the parameter NPG
. Then, the field results in the Gauß points of the selected elements are written in the transfer file. In Ex. 36, only the total strain tensor and the temperature field are transferred to HOMAT. But, note that in this .txt
file other field results could have been written by the Python script aba2hom.py
, such as the stress tensor, equivalent von Mises stress, etc. HOMAT neglects their presence and reads only the mandatory macrostrain results and in presence of thermal localization (parameter LOCATH
on keyword /PARA
of HOMAT command file) the temperature state. For each Gauß point, the script writes, the element label, then the Gauß point number followed by the field values in this point. As a 2D plane strain analysis is defined in Example 36, four strain components are written: E_{XX}, E_{YY}, E_{XY} and E_{ZZ} in the structural axis system of ABAQUS, leading to following specification:
10, 3, 0.00047, 0.01486, 0.00000, 0.00014
Remark: If the structural axis systems differ at macro and microscale, the user must specify the rotation matrix of the macrostrains in the microscale axis system. This mandatory rotation will be defined by a dedicated orientation definition via the optional parameter ORIENTATION
on the /INPUT FILES
keyword, placed there behind the MACRO
parameter, specifying the macro to micro transfer file.
Homogenization of Materials with Periodic Microstructure¶
The asymptotic homogenization method will be briefly outlined here for a thermoelastic material with a periodic microstructure subjected to volumetric eigenstrains. A more detailed and theoretical description of the method can be found in the literature ^{24}.
Static equilibrium of a heterogeneous material¶
The asymptotic homogenization method (AHM) assumes that the heterogeneous material is built at the microscale by the periodic repetition of a Representative Volume Element (RVE), At the macroscale this heterogeneous material with volume \Omega is in a static equilibrium with the body forces b_{i}, the prescribed surface tractions \overline{t}_{i} on \Gamma_{t} and the prescribed displacements \overline{u}_{i} on \Gamma_{u}. The sample is further subjected to a temperature variation \Delta T = T^{h}  T_{r}, where T_{r} is a reference temperature being set to 0°C without loss of generality.
In the macroscopic reference coordinate system \vec{z} = (z_{1}, z_{2}, z_{3}), see Fig. 25, the static equilibrium of the heterogeneous material is given by
As a pure elastic material behaviour of isotropic, orthotropic or anisotropic type is assumed for the heterogeneous material, the stress tensor \sigma^{h}_{ij} is related to the total strain \varepsilon^{h}_{kl} by following generalized Hooke law
The Hooke matrix H^{h}_{ijkl} and the thermal expansion coefficient \alpha^{h}_{kl} are both temperature dependent quantities. The volumetric eigenstrains \kappa^{h}_{kl} depend only on the position Z and are generally present in polycrystalline materials. To complete the boundary value problem, the applied surface tractions \overline{t}_{i} and displacements \overline{u}_{i} are subject to the following conditions on \Gamma_{t} and \Gamma_{u}, respectively
Periodicity and twoscale description¶
A heterogeneous material has a regular periodicity if its field variables \varphi satisfy:
with
where n_{1}, n_{2}, n_{3} are integers and \vec{Y} = (Y_{1}, Y_{2}, Y_{3}) is a constant vector, whose components describe the material periodicity (see Fig. 28).
In the homogenization theory, the periodicity \vec{Y} is small compared to the specimen's dimensions. As the function \varphi, like e.g. the Hooke matrix, varies strongly in the neighbourhood of an arbitrary point Z_{0}, (see Fig. 27), two scales are introduced in order to describe the dependency of each physical quantity:
 a macroscopic, global variable \vec{x}, which measures the variation of the mean value of function \varphi from unit cell to unit cell,
 a microscopic, periodic variable \vec{y}, which represents the strong variation of the function within each unit cell.
The scale ratio \epsilon between the micro and macroscale, defined by y=x/\epsilon, is usually a small factor, which specifies also the size ratio between the extensions of RVE and specimen. The material at the macroscale can thus be described as the combination of identical unit cells of extension (\epsilon Y_{1}, \epsilon Y_{2}, \epsilon Y_{3}) (see Fig. 28).
The asymptotic homogenization method¶
Starting point of this method is a formal asymptotic expansion of the periodic displacement and temperature fields, \vec{u}^{h} and T^{h} on the heterogeneous material with the aspect ratio \epsilon:
As the displacement and the temperature field expansions are both limited to the 1^{st} order here, the corresponding homogenization scheme is named 1^{st} order homogenization method. Expansions up to the 2^{nd} order have been introduced by Koutznetsova et al. ^{23} in order to capture microstructural size effects. This formalism is complex and leads to the introduction of an enriched macrocontinuum of Cosserat type not being discussed further here.
The Hooke matrix H^{h}_{ijkl}, the thermal expansion coefficients \alpha^{h}_{kl} and volumetric eigenstrains \kappa^{h}_{kl} are assumed to be periodic function on the unit cell and their 1^{st} order asymptotic expansions are given by:
The series expansions Eqs. \eqref{eq:asy_expansion} and (\ref{eq:asy_expansion_props}) and the differential operator
imply the following stress definitions:
where \vec{u}^{r} specifies either \vec{u}^{0} or \vec{u}^{1} and E_{kl} and e_{kl} are the macroscopic and microscopic total strains, expressed either in terms of \vec{x} or \vec{y} respectively:
These strains are introduced in the static equilibrium condition Eq. \eqref{eq:equil_condition}, leading thus to a system of differential equations, which can be sorted with respect to the scale ratio \epsilon. The differential equations related to the different powers of \epsilon are discussed in the following:

\epsilon^{2} term:
\begin{equation} \frac{\partial}{\partial y_{j}} s^{0}_{ij}\left(\vec{x},\vec{y},T^{0}\right) = 0\,. \label{eq:eps2} \end{equation}These differential equations are homogeneous. As the Hooke matrix is positive definite, e_{ij}(\vec{u}^{0}) = 0 and u^{0}(\vec{x},\vec{y}) = u(\vec{x}). Thus, the first term of the serial expansion, \vec{u}^{0}, depends only on the macroscopic variable \vec{x} and not on the microscopic one.  \epsilon^{1} term:
\begin{equation} \frac{\partial}{\partial y_{j}} S^{0}_{ij}\left(\vec{x},\vec{y},T^{0}\right) + \frac{\partial}{\partial y_{j}} s^{1}_{ij}\left(\vec{x},\vec{y},T^{0}\right)  \frac{\partial}{\partial y_{j}} p^{0}_{ij}\left(\vec{x},\vec{y},T^{0}\right) = 0\,. \label{eq:eps1} \end{equation}This differential system is linear in \vec{y}. Its solution is composed of the global solution \vec{u}^{1}(\vec{x},\vec{y})=\bar{\vec{u}}^{1}(\vec{x}) of the homogeneous system plus a special solution of type
\begin{equation} \label{eq:spec_solution} \vec{u}^{1}(\vec{x},\vec{y}) = \zeta^{rs}_{i} E^{0}_{rs}(\vec{x}) + \psi_{i}(\vec{y})\,, \end{equation}where \vec{\zeta}^{rs} = ( \zeta^{rs}_{x}, \zeta^{rs}_{y}, \zeta^{rs}_{z}), with rs = xx, yy, zz, xy, xz, yz are periodic microscopic displacement fields induced by the initial strains E^{0}_{rs}. Furthermore, \psi(\vec{y}) are microscopic displacements generated by the eigenstress tensor p^{0}_{ij} due to thermal expansion and/or volume changes of the material phases. The special solution \eqref{eq:spec_solution} is introduced into the system \eqref{eq:eps1}. Then, the initial uniform strains E^{0}_{rs}(\vec{x}) are extracted from this system, leading to:
\begin{multline} \label{eq:solu_1} E^{0}_{rs} \frac{\partial}{\partial y_{j}} \left[ H_{ijkl}\left(\vec{x},\vec{y},T^{0}\right)  H_{ijkl}\left(\vec{x},\vec{y},T^{0}\right)\left[ e_{kl}\left( \vec{\zeta}^{rs} \right) + e_{kl}\left( \vec{\psi} \right) \right] \right] + \\ + \frac{\partial}{\partial y_{j}} p^{0}_{ij}\left(\vec{x},\vec{y},T^{0}\right) = 0 \end{multline}As the initial strains and generalized eigenstrains can be specified on the unit cell separately, the modified system \eqref{eq:solu_1} can be split into:
\begin{equation} \begin{aligned} & \frac{\partial}{\partial y_{j}} H_{ijkl}\left(\vec{x},\vec{y},T^{0}\right) = H_{ijkl}\left(\vec{x},\vec{y},T^{0}\right) e_{kl}\left( \vec{\zeta}^{rs} \right)\\ & \frac{\partial}{\partial y_{j}} p^{0}_{ij}\left(\vec{x},\vec{y},T^{0}\right) = H_{ijkl}\left(\vec{x},\vec{y},T^{0}\right) e_{kl}\left( \vec{\psi} \right) \end{aligned} \label{eq:solu_2} \end{equation}A variational formulation is applied to both systems \eqref{eq:solu_2} in order to determine the unknown displacement fields \vec{\zeta}^{rs} and \vec{\psi} ^{8}. Then, the 3D RVE of the heterogeneous material is discretized by isoparametric volume elements leading to following linearized systems in the nodal values \vec{Z}^{rs} and \Psi of the unknown displacement fields \vec{\zeta}^{rs} and \vec{\psi}:
\begin{equation} \begin{aligned} & \vec{K}\left( T^{0} \right) \vec{Z}^{rs} = \vec{g}^{rs}_{\text{im,H}}\left( T^{0} \right) + \vec{g}^{rs}_{\Delta T, \text{H}}\left( T^{0} \right)\\ & \vec{K}\left( T^{0} \right) \vec{\Psi} = \vec{g}^{rs}_{\text{im,p}}\left( T^{0} \right) + \vec{g}^{rs}_{\Delta T, \text{p}}\left( T^{0} \right) \end{aligned} \label{eq:solu_3} \end{equation}where:

\vec{K}\left( T^{0} \right) is the classical elastic stiffness matrix of the RVE;

\vec{g}^{rs}_{\text{im,H}} and \vec{g}^{rs}_{\text{im,p}} are implicit surface tensions, which are induced either by the jump of the Hooke matrix or of the eigenstress at the phase/material interfaces;

\vec{g}^{rs}_{\Delta T, \text{H}} and \vec{g}^{rs}_{\Delta T, \text{p}} are implicit body forces, which are induced either by the temperature dependence of the Hooke matrix or of the eigenstress in each phase/material of the RVE.


\epsilon^{0} term: The solution of the \epsilon_{0} term leads to the definition of following effective thermoelastic properties of the heterogeneous material:

effective Hooke matrix:
\begin{equation} \label{eq:eff_hoke_mat} H^{\text{hom}}_{ijkl}\left( T^{0} \right) = \frac{1}{Y}\int_{Y}\left[ H_{ijkl}\left( y, T^{0} \right)  H_{ijrs}\left( y, T^{0} \right) e_{rs}\left( \vec{\zeta}^{kl} \right)\right]\,\text{d}y \end{equation} 
effective eigenstress:
\begin{equation} \label{eq:eff_eigenstress} p^{\text{hom}}_{ij}\left( T^{0} \right) = \frac{1}{Y}\int_{Y}\left[ p_{ij}\left( y, T^{0} \right)  H_{ijrs}\left( y, T^{0} \right) e_{rs}\left( \vec{\psi}\,, \right)\right]\,\text{d}y \end{equation}
where Y is the volume of the RVE.

The generalized effective eigenstrains \bar{\kappa}^{\text{hom}}_{ij} are obtained easily by inverting the relation:
Expressed for \Delta T^{0} = 1 Kelvin, relation \eqref{eq:gen_eff_eigenstrain} cannot simultaneously provide the effective thermal expansion coefficients \alpha^{\text{hom}}_{ij} and the effective transformation strains \kappa^{\text{hom}}_{ij}. Therefore, two effective eigenstresses are evaluated successively: one including the volumetric strain \kappa_{ij} and one considering only the thermal expansion coefficient. In addition to the purely mechanical problem (determination of the 6 periodic displacement fields \vec{\chi}^{\text{rs}}), the new 3D thermoelastic procedure will consider either 2 additional effective eigenstresses (in the presence of eigenstrain and thermal expansion) or one additional eigenstress field (only thermal expansion). This procedure is implemented in HOMAT since version 5.2 and for 3D thermoelastic homogenization the solution of either 8 or 7 linear systems.
Finally, in order to get the effective engineering properties, the predicted Hooke matrix is made symmetric and then inverted. From the effective compliance matrix the orthotropic Young (E_{x}, E_{y}, E_{z}) and shear moduli (G_{xy}, G_{xz}, G_{yz}) and the Poisson coefficients (\nu_{xy}, \nu_{xz}, \nu_{yz}) can be derived. In the special case of a quasi cubic global behavior, not only the effective Young and shear moduli and the Poisson coefficient \nu can be extracted, but also an effective Zener anisotropy coefficient, defined by
Postprocessing of macroscale results: the localization step¶
The coupling to the macroscale simulation is realized by a localization step comprising the evaluation of the microstrains and stresses on the RVE. The thermomechanical analysis using effective properties on the macroscale allows the evaluation of macrostrains and stresses in the material. In regions with high stress concentrations it is also important to accurately know the microstress and strain levels on the RVE. They are needed to predict the damage and/or failure initiation criterion on the RVE and to localize critical zones of the microstructure (see Figs. 29 and 30).
Therefore, at selected load increments the macrostrain state is extracted and transferred to the RVE. It specifies there the boundary conditions. In the framework of the asymptotic homogenization method (see Appendix Homogenization of Materials with Periodic Microstructure), the asymptotic development of the displacement field \eqref{eq:asy_expansion} is introduced in the total strain expression
The special solution \eqref{eq:spec_solution} for the displacement term \vec{u}^{1}(\vec{x},\vec{y}) allows to transform expression \eqref{eq:loc_str}, as follows
This expression specifies the total microstrains in the RVE for the given macrostrains E_{ij}^{0}(\vec{x}). Similarly, by taking into account the definitions [Eq. \eqref{eq:stresses}] of the macro and microstresses, S_{ij} and s_{ij} and the eigenstress p_{ij}^{0}\left( \vec{x},\vec{y},T \right), the microscopic stresses are given by:
Based on these local stress/strain results on the RVE, a material dependent failure criterion can be applied in order to localize failure initiation in the heterogeneous material's microstructure at a specific load level.
Periodic boundary conditions via a masterslave projection algorithm¶
The application of the periodic boundary conditions on the RVE via the penalty method works but leads to a poorly conditioned system matrix with high differences between the stiffness terms and induces locally additional artificial stiffness. In order to more accurately predict effective mechanical properties, a masterslave projection algorithm, detailed in Zhang and Oskay ^{25}, has been implemented in HOMAT version 6.0. Complex RVE's, like the crosshatched PP bilamella mesh are often at best close to but not strictly periodic. In strictly periodic microstructures, the periodicity is simply enforced by assigning masterslave pairs on opposite RVE surfaces and eliminate the slave nodes from the global systems \eqref{eq:solu_3}, as they have the same displacements at the corresponding master nodes. For nonmatching meshes on opposite periodic surfaces, we specify the surface with the finer mesh density to be the slave surface and the other one to be the master one. Then, each node, belonging to the slave surface element is projected on the master surface (see Fig. 31). Its specific projection point in a masterelement is identified, allowing its interpolation in this specific element. Performing this nodal projection for all nodes on a slave surface, yields the following relationship between slave and master displacement vectors
where \overline{Z}_{s}^{rs,e} and \overline{Z}_{m}^{rs,e} denote the vectors of slave and master nodes respectively and \overline{C} is the rectangular projection matrix build by the local interpolation coefficients. Expression \eqref{eq:mast_slave} is then used to eliminate the slave nodes per surface element of the slave surface. The global matrix systems \eqref{eq:solu_3} are thus reduced and contain only bulk and surface master nodes.

Georg J. Schmitz and Ulrich Prahl. Integrative Computational Materials Engineering: Concepts and Applications of a Modular Simulation Platform. John Wiley & Sons, July 2012. ISBN 9783527646111. ↩↩

Gottfried Laschet, Hamed Nokhostin, Simon Koch, Marc Meunier, and Christian Hopmann. Prediction of effective elastic properties of a polypropylene component by an enhanced multiscale simulation of the injection molding process. Mechanics of Materials, 140:103225, January 2020. doi:10.1016/j.mechmat.2019.103225. ↩↩↩

Gottfried Laschet, Marcel Spekowius, Roberto Spina, and Christian Hopmann. Multiscale simulation to predict microstructure dependent effective elastic properties of an injection molded polypropylene component. Mechanics of Materials, 105:123–137, February 2017. doi:10.1016/j.mechmat.2016.10.009. ↩↩↩↩

Marcel Spekowius, Roberto Spina, and Christian Hopmann. Mesoscale simulation of the solidification process in injection moulded parts. Journal of Polymer Engineering, 366:563–573, August 2016. doi:10.1515/polyeng20140223. ↩

K. Bobzin, N. Bagcivan, D. Parkot, T. Kashko, G. Laschet, and J. Scheele. Influence of the Definition of the Representative Volume Element on Effective Thermoelastic Properties of Thermal Barrier Coatings with Random Microstructure. Journal of Thermal Spray Technology, 185:988, July 2009. doi:10.1007/s1166600993510. ↩↩↩↩

Dassault Systèmes Simulia Corp. Abaqus Finite Element Analysis, User manual. Providence, RI, USA, 2019. ↩↩↩

Gottfried Laschet. Homogenization of the thermal properties of transpiration cooled multilayer plates. Computer Methods in Applied Mechanics and Engineering, 19141:4535–4554, September 2002. doi:10.1016/S0045782502003195. ↩

Gottfried Laschet. Forchheimer law derived by homogenization of gas flow in turbomachines. Journal of Computational and Applied Mathematics, 2152:467–476, June 2008. doi:10.1016/j.cam.2006.03.049. ↩↩↩↩

Jian Li and Yinnian He. A stabilized finite element method based on two local Gauss integrations for the Stokes equations. Journal of Computational and Applied Mathematics, 2141:58–65, April 2008. doi:10.1016/j.cam.2007.02.015. ↩↩

O C Zienkiewicz and R. L. Taylor. The Finite Element Method 4/e. Mc Graw Hill, Singapore, 1989. ISBN 9780070991934. ↩↩↩

Gottfried Laschet, Jörg Sauerhering, Oliver Reutter, Thomas Fend, and Josef Scheele. Effective permeability and thermal conductivity of opencell metallic foams via homogenization on a microstructure model. Computational Materials Science, 453:597–603, May 2009. doi:10.1016/j.commatsci.2008.06.023. ↩↩

G. Laschet, P. Fayek, T. Henke, H. Quade, and U. Prahl. Derivation of anisotropic flow curves of ferrite–pearlite pipeline steel via a twolevel homogenisation scheme. Materials Science and Engineering: A, 566:143–156, March 2013. doi:10.1016/j.msea.2012.12.064. ↩

S. NematNasser and M. Hori. Micromechanics: Overall Properties of Heterogeneous Materials. Elsevier, October 2013. ISBN 9781483291512. ↩

JeanMarie Berthelot. Composite Materials: Mechanical Behavior and Structural Analysis. Springer, New York, 1999 edition edition, December 1998. ISBN 9780387984261. ↩↩

D. N. Arnold, F. Brezzi, and M. Fortin. A stable finite element for the stokes equations. CALCOLO, 214:337–344, December 1984. doi:10.1007/BF02576171. ↩

Gottfried Laschet and Stephan Rex. Effective Thermoelastic Properties of Flat and Curved Transpiration Cooled Multilayer Plates via Homogenization. In ASME Turbo Expo 2008: Power for Land, Sea, and Air, 1757–1768. American Society of Mechanical Engineers Digital Collection, August 2009. doi:10.1115/GT200850590. ↩

Claude Pommerell. Solution of Large Unsymmetric Systems of Linear Equations. PhD thesis, ETH Zurich, 1992. doi:10.3929/ETHZA000669614. ↩

Jan Dvorak. A Reliable Numerical Method for Computing Homogenized Coefficients. Mathematical Institute, Danish Technical University, Lyngby, Denmark, 1995. ↩

JeanBaptiste Clavaud, Alexis Maineult, Maria Zamora, Patrick Rasolofosaon, and Camille Schlitter. Permeability anisotropy and its relations with porous medium structure. Journal of Geophysical Research: Solid Earth, 2008. _eprint: https://agupubs.onlinelibrary.wiley.com/doi/pdf/10.1029/2007JB005004. doi:10.1029/2007JB005004. ↩

A. Tamayol and M. Bahrami. Parallel Flow Through Ordered Fibers: An Analytical Approach. Journal of Fluids Engineering, November 2010. doi:10.1115/1.4002169. ↩

Ali Tamayol and Majid Bahrami. Transverse permeability of fibrous porous media. Physical Review E, 834:046314, April 2011. doi:10.1103/PhysRevE.83.046314. ↩

Gottfried Laschet, Mohit Shukla, Thomas Henke, Patrick Fayek, Markus Bambach, and Ulrich Prahl. Impact of the Microstructure on the U–O Forming Simulations of a Ferrite–Pearlite Pipeline Tube. steel research international, 856:1083–1098, 2014. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/srin.201300247. doi:10.1002/srin.201300247. ↩

V. Kouznetsova, M. G. D. Geers, and W. a. M. Brekelmans. Multiscale constitutive modelling of heterogeneous materials with a gradientenhanced computational homogenization scheme. International Journal for Numerical Methods in Engineering, 548:1235–1260, 2002. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/nme.541. doi:10.1002/nme.541. ↩↩

JoséMiranda Guedes and Noboru Kikuchi. Preprocessing and postprocessing for materials based on the homogenization method with adaptive finite element methods. Computer Methods in Applied Mechanics and Engineering, 832:143–198, October 1990. doi:10.1016/004578259090148F. ↩

Xiang Zhang and Caglar Oskay. Eigenstrain based reduced order homogenization for polycrystalline materials. Computer Methods in Applied Mechanics and Engineering, 297:408–436, December 2015. doi:10.1016/j.cma.2015.09.006. ↩

The material properties are tabulated as a function of a scalar T, which is typically interpreted as temperature. However, you are free to change this interpretation according to your needs, e.g., if you prefer a description based on the fraction solid or concentration, you can provide material data as a function of these quantities and supply a fraction solid or concentration wherever HOMAT requires a temperature. ↩↩↩↩↩↩↩↩↩↩

Note that a material is not detected as cubic, if the cubic axes do not coincide with the RVE axes. ↩