40 #ifdef _WIN32 //_MSC_VER and __MINGW32__ included 55 umatobj(NULL), umat(NULL),
56 mStressInterpretation(0),
57 mUseNumericalTangent(false),
65 FreeLibrary( ( HMODULE ) this->
umatobj );
78 std :: string umatname;
81 if ( result !=
IRRT_OK )
return result;
88 strncpy(this->
cmname, umatname.c_str(), 80);
95 DWORD dlresult = GetLastError();
100 * ( FARPROC * ) ( & this->
umat ) = GetProcAddress( ( HMODULE ) this->
umatobj,
"umat_" );
103 DWORD dlresult = GetLastError();
104 OOFEM_ERROR(
"Couldn't load symbol umat,\nerror code: %d\n", dlresult);
108 this->umatobj = dlopen(
filename.c_str(), RTLD_NOW);
109 if ( !this->umatobj ) {
113 * (
void ** ) ( & this->
umat ) = dlsym(this->umatobj,
"umat_");
115 char *dlresult = dlerror();
117 OOFEM_ERROR(
"couldn't load symbol umat,\ndlerror: %s\n", dlresult);
169 for (
int i = 1; i <= strain.
giveSize(); ++i ) {
174 stressh.
times(1.0 / h);
175 En.setColumn(stressh, i);
181 printf(
"Tangent = ");
217 for (
int i = 1; i <= 9; ++i ) {
222 stressh.
times(1.0 / h);
268 int ntens = ndi + nshr;
280 double temp = 0.0, dtemp = 0.0;
285 double time [ 2 ] = {
294 double sse, spd, scd;
355 OOFEM_LOG_DEBUG(
"AbaqusUserMaterial :: giveRealStressVector_3d - Calling subroutine");
408 OOFEM_LOG_DEBUG(
"AbaqusUserMaterial :: giveRealStressVector - Calling subroutine was successful");
425 int ntens = ndi + nshr;
447 double temp = 0.0, dtemp = 0.0;
452 double time [ 2 ] = {
461 double sse, spd, scd;
522 OOFEM_LOG_DEBUG(
"AbaqusUserMaterial :: giveRealStressVector - Calling subroutine");
616 OOFEM_LOG_DEBUG(
"AbaqusUserMaterial :: giveRealStressVector_3d - Calling subroutine was successful");
623 if ( type == IST_Undefined || type == IST_AbaqusStateVector ) {
627 }
else if ( type == IST_StressTensor ) {
641 tempStateVector = stateVector;
646 numState(numState), stateVector(numState), tempStateVector(numState), hasTangentFlag(false)
663 fprintf(File,
" stateVector ");
665 for (
auto &var : state ) {
666 fprintf( File,
" % .4e", var );
double giveDeterminant() const
Returns the trace (sum of diagonal components) of the receiver.
virtual void updateYourself(TimeStep *tStep)
Update equilibrium history variables according to temp-variables.
InternalStateType
Type representing the physical meaning of element or constitutive model internal variable.
void subtract(const FloatArray &src)
Subtracts array src to receiver.
#define _IFT_AbaqusUserMaterial_initialStress
void letTempStrainVectorBe(const FloatArray &v)
Assigns tempStrainVector to given vector v.
virtual MaterialStatus * giveStatus(GaussPoint *gp) const
Returns material status of receiver in given integration point.
virtual void giveRealStressVector_3d(FloatArray &answer, GaussPoint *gp, const FloatArray &reducedStrain, TimeStep *tStep)
Default implementation relies on giveRealStressVector for second Piola-Kirchoff stress.
void letTempFVectorBe(const FloatArray &v)
Assigns tempFVector to given vector v.
double mPerturbation
Size of perturbation if numerical tangent is used.
FloatArray & letTempStateVectorBe(FloatArray &s)
AbaqusUserMaterial(int n, Domain *d)
Constructor.
FloatArray tempStateVector
Temporary state vector.
virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in Reduced form.
double & at(int i)
Coefficient access function.
int numState
Size of the state vector.
This class implements a structural material status information.
char cmname[80]
Name for material routine.
virtual void printOutputAt(FILE *file, TimeStep *tStep)
Print receiver's output to given stream.
virtual void give3dMaterialStiffnessMatrix_dPdF(FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
void beVectorForm(const FloatMatrix &aMatrix)
Reciever will be a vector with 9 components formed from a 3x3 matrix.
virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in Reduced form.
double giveTargetTime()
Returns target time.
virtual void giveInputRecord(DynamicInputRecord &input)
Setups the input record string of receiver.
Element * giveElement()
Returns corresponding element to receiver.
AbaqusUserMaterialStatus(int n, Domain *d, GaussPoint *gp, int numState)
Constructor.
const double * givePointer() const
Exposes the internal values of the matrix.
#define _IFT_AbaqusUserMaterial_numericalTangentPerturbation
virtual void printOutputAt(FILE *file, TimeStep *tStep)
Print receiver's output to given stream.
#define OOFEM_LOG_DEBUG(...)
MatResponseMode
Describes the character of characteristic material matrix.
virtual void give3dMaterialStiffnessMatrix(FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
Computes full 3d material stiffness matrix at given integration point, time, respecting load history ...
#define _IFT_AbaqusUserMaterial_numState
void changeComponentOrder()
Swaps the fourth and sixth index in the array.
void beDifferenceOf(const FloatArray &a, const FloatArray &b)
Sets receiver to be a - b.
double giveTimeIncrement()
Returns solution step associated time increment.
static int n
Gausspoint counter.
void beMatrixForm(const FloatArray &aArray)
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
FloatArray stateVector
General state vector.
bool mUseNumericalTangent
Flag to determine if numerical tangent should be used.
#define _IFT_AbaqusUserMaterial_numericalTangent
const FloatArray & giveStateVector() const
void times(double f)
Multiplies receiver by factor f.
virtual void giveFirstPKStressVector_3d(FloatArray &answer, GaussPoint *gp, const FloatArray &reducedF, TimeStep *tStep)
Default implementation relies on giveRealStressVector for second Piola-Kirchoff stress.
std::string filename
Name of the file that contains the umat function.
double at(int i, int j) const
Coefficient access function.
#define _IFT_AbaqusUserMaterial_name
virtual ~AbaqusUserMaterial()
Destructor.
virtual IRResultType initializeFrom(InputRecord *ir)
Reads the following values;.
Abstract base class representing a material status information.
virtual void initTempStatus()
Initializes the temporary internal variables, describing the current state according to previously re...
Class representing vector of real numbers.
int mStressInterpretation
Flag to determine how the stress and Jacobian are interpreted.
FloatArray initialStress
Initial stress.
FloatArray strainVector
Equilibrated strain vector in reduced form.
Implementation of matrix containing floating point numbers.
IRResultType
Type defining the return values of InputRecord reading operations.
FloatArray properties
Material properties.
void changeComponentOrder()
Swaps the indices in the 6x6 matrix such that it converts between OOFEM's and Abaqus' way of writing ...
void letTempStressVectorBe(const FloatArray &v)
Assigns tempStressVector to given vector v.
const FloatArray & giveStressVector() const
Returns the const pointer to receiver's stress vector.
IntegrationPointStatus * giveMaterialStatus()
Returns reference to associated material status (NULL if not defined).
void * umatobj
Dynamically loaded umat.
const FloatArray & giveFVector() const
Returns the const pointer to receiver's deformation gradient vector.
int giveMetaStepNumber()
Returns receiver's meta step number.
static void giveReducedMatrixForm(FloatMatrix &answer, const FloatMatrix &full, MaterialMode matMode)
Converts the full symmetric Voigt matrix (4th order tensor) to reduced form.
void zero()
Zeroes all coefficients of receiver.
void beTProductOf(const FloatMatrix &a, const FloatMatrix &b)
Assigns to the receiver product of .
const FloatArray & giveTempFVector() const
Returns the const pointer to receiver's temporary deformation gradient vector.
void setColumn(const FloatArray &src, int c)
Sets the values of the matrix in specified column.
void times(double s)
Multiplies receiver with scalar.
const double * givePointer() const
Gives the pointer to the raw data, breaking encapsulation.
virtual void updateYourself(TimeStep *tStep)
Update equilibrium history variables according to temp-variables.
void printYourself() const
Prints matrix to stdout. Useful for debugging.
void beSymVectorForm(const FloatMatrix &aMatrix)
Reciever will be a vector with 6 components formed from a 3x3 matrix.
Abstract base class for all "structural" constitutive models.
void beProductTOf(const FloatMatrix &a, const FloatMatrix &b)
Assigns to the receiver product of .
virtual void initTempStatus()
Initializes the temporary internal variables, describing the current state according to previously re...
Domain * giveDomain() const
void beUnitMatrix()
Sets receiver to unity matrix.
virtual void giveInputRecord(DynamicInputRecord &input)
Setups the input record string of receiver.
REGISTER_Material(DummyMaterial)
void beSymVectorFormOfStrain(const FloatMatrix &aMatrix)
void beProductOf(const FloatMatrix &a, const FloatMatrix &b)
Assigns to the receiver product of .
int giveSize() const
Returns the size of receiver.
the oofem namespace is to define a context or scope in which all oofem names are defined.
#define _IFT_AbaqusUserMaterial_properties
void(* umat)(double *stress, double *statev, double *ddsdde, double *sse, double *spd, double *scd, double *rpl, double *ddsddt, double *drplde, double *drpldt, double *stran, double *dstran, double time[2], double *dtime, double *temp, double *dtemp, double predef[1], double dpred[1], char cmname[80], int *ndi, int *nshr, int *ntens, int *nstatv, double *props, int *nprops, double coords[3], double *drot, double *pnewdt, double *celent, double *dfgrd0, double *dfgrd1, int *noel, int *npt, int *layer, int *kspt, int *kstep, int *kinc)
Pointer to the dynamically loaded umat-function (translated to C)
void beInverseOf(const FloatMatrix &src)
Modifies receiver to become inverse of given parameter.
const FloatArray & giveStrainVector() const
Returns the const pointer to receiver's strain vector.
virtual MaterialStatus * CreateStatus(GaussPoint *gp) const
Creates new copy of associated status and inserts it into given integration point.
Class representing integration point in finite element program.
#define _IFT_AbaqusUserMaterial_userMaterial
Class representing solution step.
virtual void givePlaneStrainStiffMtrx_dPdF(FloatMatrix &answer, MatResponseMode mmode, GaussPoint *gp, TimeStep *tStep)
virtual int computeGlobalCoordinates(FloatArray &answer, const FloatArray &lcoords)
Computes the global coordinates from given element's local coordinates.
void add(const FloatArray &src)
Adds array src to receiver.
void letTempTangentBe(FloatMatrix t)
const FloatArray & giveNaturalCoordinates()
Returns coordinate array of receiver.
void letTempPVectorBe(const FloatArray &v)
Assigns tempPVector to given vector v.
void resize(int s)
Resizes receiver towards requested size.
const FloatMatrix & giveTempTangent()