v0.13.1
Classes | Functions
SmallStrainPlasticity.hpp File Reference

Operators and data structures for small strain plasticity. More...

Go to the source code of this file.

Classes

struct  SmallStrainPlasticity
 Small strain finite element implementation. More...
 
struct  SmallStrainPlasticity::MyVolumeFE
 definition of volume element More...
 
struct  SmallStrainPlasticity::CommonData
 common data used by volume elements More...
 
struct  SmallStrainPlasticity::OpGetDataAtGaussPts
 
struct  SmallStrainPlasticity::OpGetCommonDataAtGaussPts
 
struct  SmallStrainPlasticity::OpUpdate
 
struct  SmallStrainPlasticity::MakeB
 
struct  SmallStrainPlasticity::OpAssembleRhs
 
struct  SmallStrainPlasticity::OpAssembleLhs
 
struct  SmallStrainPlasticity::ClosestPointProjection
 Closest point projection algorithm. More...
 
struct  SmallStrainPlasticity::OpCalculateStress
 

Functions

PetscErrorCode SmallStrainPlasticityfRes (SNES snes, Vec chi, Vec r, void *ctx)
 
PetscErrorCode SmallStrainPlasticityfJac (SNES, Vec, Mat, Mat, void *ctx)
 

Detailed Description

Operators and data structures for small strain plasticity.

Definition in file SmallStrainPlasticity.hpp.

Function Documentation

◆ SmallStrainPlasticityfJac()

PetscErrorCode SmallStrainPlasticityfJac ( SNES  snes,
Vec  chi,
Mat  A,
Mat  ,
void *  ctx 
)
Examples
SmallStrainPlasticity.hpp.

Definition at line 554 of file SmallStrainPlasticity.cpp.

1192 {
1193 PetscFunctionBegin;
1196
1197 {
1198 #if PETSC_VERSION_GE(3,6,0)
1199 const double *array;
1200 ierr = VecGetArrayRead(chi,&array); CHKERRQ(ierr);
1201 #else
1202 double *array;
1203 ierr = VecGetArray(chi,&array); CHKERRQ(ierr);
1204 #endif
1205 for(int ii = 0;ii<6;ii++) {
1206 cp->plasticStrain[ii] = array[ii];
1207 }
1208 for(unsigned int ii = 0;ii<cp->internalVariables.size();ii++) {
1209 cp->internalVariables[ii] = array[6+ii];
1210 }
1211 cp->deltaGamma = array[6+cp->internalVariables.size()];
1212 #if PETSC_VERSION_GE(3,6,0)
1213 ierr = VecRestoreArrayRead(chi,&array); CHKERRQ(ierr);
1214 #else
1215 ierr = VecRestoreArray(chi,&array); CHKERRQ(ierr);
1216 #endif
1217 }
1218 ierr = cp->evaluatePotentials(); CHKERRQ(ierr);
1219 ierr = cp->cAlculateA(); CHKERRQ(ierr);
1220 PetscFunctionReturn(0);
1221}
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76

◆ SmallStrainPlasticityfRes()

PetscErrorCode SmallStrainPlasticityfRes ( SNES  snes,
Vec  chi,
Vec  r,
void *  ctx 
)
Examples
SmallStrainPlasticity.hpp.

Definition at line 553 of file SmallStrainPlasticity.cpp.

1139 {
1140 PetscFunctionBegin;
1143
1144 //ierr = VecView(chi,PETSC_VIEWER_STDOUT_WORLD); CHKERRQ(ierr);
1145 {
1146 #if PETSC_VERSION_GE(3,6,0)
1147 const double *array;
1148 ierr = VecGetArrayRead(chi,&array); CHKERRQ(ierr);
1149 #else
1150 const double *array;
1151 ierr = VecGetArray(chi,&array); CHKERRQ(ierr);
1152 #endif
1153 for(int ii = 0;ii<6;ii++) {
1154 cp->plasticStrain[ii] = array[ii];
1155 }
1156 for(unsigned int ii = 0;ii<cp->internalVariables.size();ii++) {
1157 cp->internalVariables[ii] = array[6+ii];
1158 }
1159 cp->deltaGamma = array[6+cp->internalVariables.size()];
1160 #if PETSC_VERSION_GE(3,6,0)
1161 ierr = VecRestoreArrayRead(chi,&array); CHKERRQ(ierr);
1162 #else
1163 ierr = VecRestoreArray(chi,&array); CHKERRQ(ierr);
1164 #endif
1165 }
1166 ierr = cp->evaluatePotentials(); CHKERRQ(ierr);
1167 ierr = cp->cAlculateR(r); CHKERRQ(ierr);
1168 // ierr = VecView(r,PETSC_VIEWER_STDOUT_WORLD); CHKERRQ(ierr);
1169 // cout << "sTress " << cp->sTress << endl;
1170 // cout << "internalVariables " << cp->internalVariables << endl;
1171 // cout << "C " << cp->C << endl;
1172 // cout << "D " << cp->D << endl;
1173 // cout << "y " << cp->y << endl;
1174 // cout << "h " << cp->h << endl;
1175 // cout << "partialYSigma " << cp->partialYSigma << endl;
1176 // cout << "partialYFlux " << cp->partialYFlux << endl;
1177 // cout << "partialHSigma " << cp->partialHSigma << endl;
1178 // cout << "partialHFlux " << cp->partialHFlux << endl;
1179 // cout << "partial2HSigma " << cp->partial2HSigma << endl;
1180 // cout << "partial2HFlux " << cp->partial2HFlux << endl;
1181 // cout << "partial2HSigmaFlux " << cp->partial2HSigmaFlux << endl;
1182 // cout << "C " << cp->C << endl;
1183 // cout << "D " << cp->D << endl;
1184 // cout << "partialWStrainPlasticStrain " << cp->partialWStrainPlasticStrain << endl;
1185 // cout << "plasticStrain " << cp->plasticStrain << endl;
1186 // cout << "internalVariables " << cp->internalVariables << endl;
1187 // cout << "delta plasticStrain " << cp->plasticStrain-cp->plasticStrain0 << endl;
1188 // cout << "delta internalVariables " << cp->internalVariables-cp->internalVariables << endl;
1189 PetscFunctionReturn(0);
1190}
const double r
rate factor