v0.13.1
Public Member Functions | Public Attributes | Friends | List of all members
SmallStrainPlasticity::ClosestPointProjection Struct Reference

Closest point projection algorithm. More...

#include <users_modules/small_strain_plasticity/src/SmallStrainPlasticity.hpp>

Inheritance diagram for SmallStrainPlasticity::ClosestPointProjection:
[legend]
Collaboration diagram for SmallStrainPlasticity::ClosestPointProjection:
[legend]

Public Member Functions

 ClosestPointProjection ()
 
virtual PetscErrorCode setActiveVariablesW ()
 
virtual PetscErrorCode setActiveVariablesYH ()
 
virtual PetscErrorCode rEcordW ()
 
virtual PetscErrorCode rEcordY ()
 
virtual PetscErrorCode rEcordH ()
 
virtual PetscErrorCode pLayW ()
 
virtual PetscErrorCode pLayW_NoHessian ()
 
virtual PetscErrorCode pLayY ()
 
virtual PetscErrorCode pLayY_NoGradient ()
 
virtual PetscErrorCode pLayH ()
 
virtual PetscErrorCode createMatAVecR ()
 
virtual PetscErrorCode destroyMatAVecR ()
 
virtual PetscErrorCode evaluatePotentials ()
 
virtual PetscErrorCode cAlculateR (Vec R)
 
virtual PetscErrorCode cAlculateA ()
 
PetscErrorCode snesCreate ()
 
PetscErrorCode snesDestroy ()
 
PetscErrorCode solveColasetProjection ()
 
PetscErrorCode consistentTangent ()
 
virtual PetscErrorCode freeHemholtzEnergy ()
 
virtual PetscErrorCode yieldFunction ()
 
virtual PetscErrorCode flowPotential ()
 

Public Attributes

VectorDouble sTrain
 
VectorDouble internalVariables
 
VectorDouble internalVariables0
 
VectorDouble plasticStrain0
 
VectorDouble plasticStrain
 
double deltaGamma
 
double tOl
 
int gG
 
VectorAdaptor sTress
 
VectorAdaptor internalFluxes
 
ublas::symmetric_matrix< double, ublas::lower > C
 
ublas::symmetric_matrix< double, ublas::lower > D
 
MatrixDouble partialWStrainPlasticStrain
 
VectorAdaptor partialYSigma
 
VectorAdaptor partialYFlux
 
VectorAdaptor partialHSigma
 
VectorAdaptor partialHFlux
 
ublas::symmetric_matrix< double, ublas::lower > partial2HSigma
 
ublas::symmetric_matrix< double, ublas::lower > partial2HFlux
 
MatrixDouble partial2HSigmaFlux
 
vector< int > tAgs
 
VectorDouble activeVariablesW
 
VectorDouble activeVariablesYH
 
double w
 
double y
 
double h
 
ublas::vector< adoublea_sTrain
 
ublas::vector< adoublea_plasticStrain
 
ublas::vector< adoublea_internalVariables
 
ublas::vector< adoublea_sTress
 
ublas::vector< adoublea_internalFluxes
 
adouble a_w
 
adouble a_y
 
adouble a_h
 
VectorDouble gradientW
 
MatrixDouble hessianW
 
VectorDouble gradientY
 
VectorDouble gradientH
 
MatrixDouble hessianH
 
Mat A
 
ublas::matrix< double, ublas::column_major > dataA
 
Vec R
 
Vec Chi
 
Vec dChi
 
SNES sNes
 
KSP kSp
 
PC pC
 
SNESLineSearch lineSearch
 
MatrixDouble Ep
 
MatrixDouble Lp
 
MatrixDouble Cp
 
MatrixDouble Cep
 
VectorDouble Yp
 

Friends

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

Detailed Description

Closest point projection algorithm.

Examples
SmallStrainPlasticityMaterialModels.hpp.

Definition at line 482 of file SmallStrainPlasticity.hpp.

Constructor & Destructor Documentation

◆ ClosestPointProjection()

SmallStrainPlasticity::ClosestPointProjection::ClosestPointProjection ( )

Member Function Documentation

◆ cAlculateA()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::cAlculateA ( )
virtual
Examples
SmallStrainPlasticity.hpp.

Definition at line 963 of file SmallStrainPlasticity.cpp.

963 {
964 PetscFunctionBegin;
965 // matrix A
966 ierr = MatZeroEntries(A); CHKERRQ(ierr);
967 // row 0
970 for(int ii = 0;ii<6;ii++) {
971 for(int jj = 0;jj<6;jj++) {
972 dataA(ii,jj) = deltaGamma*a00(ii,jj);
973 }
974 for(unsigned int jj = 0;jj<internalVariables.size();jj++) {
975 dataA(ii,6+jj) = deltaGamma*a01(ii,jj);
976 }
977 dataA(ii,ii) -= 1;
978 }
979 for(int ii = 0;ii<6;ii++) {
980 dataA(ii,6+internalVariables.size()) = partialHSigma[ii];
981 }
982 // row 1
983 MatrixDouble a11 = prod(partial2HFlux,D);
985 for(unsigned int ii = 0;ii<internalVariables.size();ii++) {
986 for(int jj = 0;jj<6;jj++) {
987 dataA(6+ii,jj) += deltaGamma*a10(ii,jj);
988 }
989 for(unsigned int jj = 0;jj<internalVariables.size();jj++) {
990 dataA(6+ii,6+jj) += deltaGamma*a11(ii,jj);
991 }
992 dataA(6+ii,6+ii) -= 1;
993 }
994 for(unsigned int ii = 0;ii<internalVariables.size();ii++) {
995 dataA(6+ii,6+internalVariables.size()) = partialHFlux[ii];
996 }
997 // row 3
998 VectorDouble partial_y_sigma_c = prod(trans(partialWStrainPlasticStrain),partialYSigma);
999 VectorDouble partial_y_flux_d = prod(trans(D),partialYFlux);
1000 for(unsigned int dd = 0;dd<partial_y_sigma_c.size();dd++) {
1001 dataA(6+internalVariables.size(),dd) = partial_y_sigma_c[dd];
1002 }
1003 for(unsigned int dd = 0;dd<partial_y_flux_d.size();dd++) {
1004 dataA(6+internalVariables.size(),6+dd) = partial_y_flux_d[dd];
1005 }
1006 dataA(6+internalVariables.size(),6+internalVariables.size()) = 0;
1007 ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
1008 ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
1009 PetscFunctionReturn(0);
1010}
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:77
UBlasVector< double > VectorDouble
Definition: Types.hpp:68
ublas::matrix< double, ublas::column_major > dataA
ublas::symmetric_matrix< double, ublas::lower > partial2HSigma
ublas::symmetric_matrix< double, ublas::lower > partial2HFlux
ublas::symmetric_matrix< double, ublas::lower > D

◆ cAlculateR()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::cAlculateR ( Vec  R)
virtual
Examples
SmallStrainPlasticity.hpp.

Definition at line 947 of file SmallStrainPlasticity.cpp.

947 {
948 PetscFunctionBegin;
949 // Residual
950 double *array;
951 ierr = VecGetArray(R,&array); CHKERRQ(ierr);
952 for(int ii = 0;ii<6;ii++) {
953 array[ii] = -plasticStrain[ii] + plasticStrain0[ii] + deltaGamma*partialHSigma[ii];
954 }
955 for(unsigned int ii = 0;ii<internalVariables.size();ii++) {
956 array[6+ii] = -internalVariables[ii] + internalVariables0[ii] + deltaGamma*partialHFlux[ii];
957 }
958 array[6+internalVariables.size()] = y;
959 ierr = VecRestoreArray(R,&array); CHKERRQ(ierr);
960 PetscFunctionReturn(0);
961}

◆ consistentTangent()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::consistentTangent ( )
Examples
SmallStrainPlasticity.hpp.

Definition at line 1080 of file SmallStrainPlasticity.cpp.

1080 {
1081 PetscFunctionBegin;
1082 ierr = evaluatePotentials(); CHKERRQ(ierr);
1083 ierr = cAlculateA(); CHKERRQ(ierr);
1084 Ep.resize(6,6,false);
1085 Lp.resize(internalVariables.size(),6,false);
1086 Yp.resize(6,false);
1087 MatrixDouble ep_row;
1088 ep_row = deltaGamma*prod(partial2HSigma,C);
1089 MatrixDouble alpha_row;
1090 alpha_row = deltaGamma*prod(trans(partial2HSigmaFlux),C);
1091 VectorDouble y_row;
1092 y_row = prod(trans(C),partialYSigma);
1093 const int n = 6+internalVariables.size();
1094 for(int dd = 0;dd<6;dd++) {
1095 ierr = VecZeroEntries(R); CHKERRQ(ierr);
1096 double *array;
1097 ierr = VecGetArray(R,&array); CHKERRQ(ierr);
1098 {
1099 for(int ii = 0;ii<6;ii++) {
1100 array[ii] = ep_row(ii,dd);
1101 }
1102 for(unsigned int ii = 0;ii<internalVariables.size();ii++) {
1103 array[6+ii] = alpha_row(ii,dd);
1104 }
1105 array[n] = y_row[dd];
1106 }
1107 ierr = VecRestoreArray(R,&array); CHKERRQ(ierr);
1108 ierr = VecAssemblyBegin(R); CHKERRQ(ierr);
1109 ierr = VecAssemblyEnd(R); CHKERRQ(ierr);
1110 ierr = KSPSolve(kSp,R,dChi); CHKERRQ(ierr);
1111 ierr = VecGetArray(dChi,&array); CHKERRQ(ierr);
1112 {
1113 for(int ii = 0;ii<6;ii++) {
1114 Ep(ii,dd) = array[ii];
1115 }
1116 for(unsigned int ii = 0;ii<internalVariables.size();ii++) {
1117 Lp(ii,dd) = array[6+ii];
1118 }
1119 Yp[dd] = array[n];
1120 }
1121 ierr = VecRestoreArray(dChi,&array); CHKERRQ(ierr);
1122 }
1123 Cp.resize(6,6,false);
1124 noalias(Cp) = prod(C,Ep);
1125 Cep.resize(6,6,false);
1126 noalias(Cep) = C+Cp;
1127 // cout << endl;
1128 // cout << "A " << dataA << endl;
1129 // cout << "C " << C << endl;
1130 // cout << "Ep " << Ep << endl;
1131 // // cout << "Lp " << Lp << endl;
1132 // // cout << "Yp " << Yp << endl;
1133 // cout << "Cp " << Cp << endl;
1134 // cout << "Cep " << Cep << endl;
1135 PetscFunctionReturn(0);
1136}
FTensor::Index< 'n', SPACE_DIM > n
ublas::symmetric_matrix< double, ublas::lower > C

◆ createMatAVecR()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::createMatAVecR ( )
virtual
Examples
SmallStrainPlasticity.hpp.

Definition at line 909 of file SmallStrainPlasticity.cpp.

909 {
910 PetscFunctionBegin;
911 PetscInt n;
912 n = 1+6+internalVariables.size();
913 dataA.resize(n,n,false);
914 ierr = MatCreateSeqDense(PETSC_COMM_SELF,n,n,&dataA(0,0),&A); CHKERRQ(ierr);
915 ierr = VecCreateSeq(PETSC_COMM_SELF,n,&R); CHKERRQ(ierr);
916 ierr = VecCreateSeq(PETSC_COMM_SELF,n,&Chi); CHKERRQ(ierr);
917 ierr = VecCreateSeq(PETSC_COMM_SELF,n,&dChi); CHKERRQ(ierr);
918 PetscFunctionReturn(0);
919}

◆ destroyMatAVecR()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::destroyMatAVecR ( )
virtual
Examples
SmallStrainPlasticity.hpp.

Definition at line 921 of file SmallStrainPlasticity.cpp.

921 {
922 PetscFunctionBegin;
923 ierr = MatDestroy(&A); CHKERRQ(ierr);
924 ierr = VecDestroy(&R); CHKERRQ(ierr);
925 ierr = VecDestroy(&Chi); CHKERRQ(ierr);
926 ierr = VecDestroy(&dChi); CHKERRQ(ierr);
927 PetscFunctionReturn(0);
928}

◆ evaluatePotentials()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::evaluatePotentials ( )
virtual
Examples
SmallStrainPlasticity.hpp.

Definition at line 930 of file SmallStrainPlasticity.cpp.

930 {
931 PetscFunctionBegin;
932 if(gG == 0) {
933 ierr = rEcordW(); CHKERRQ(ierr);
934 }
935 ierr = setActiveVariablesW(); CHKERRQ(ierr);
936 ierr = pLayW(); CHKERRQ(ierr);
937 if(gG == 0) {
938 ierr = rEcordY(); CHKERRQ(ierr);
939 ierr = rEcordH(); CHKERRQ(ierr);
940 }
941 ierr = setActiveVariablesYH(); CHKERRQ(ierr);
942 ierr = pLayY(); CHKERRQ(ierr);
943 ierr = pLayH(); CHKERRQ(ierr);
944 PetscFunctionReturn(0);
945}

◆ flowPotential()

virtual PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::flowPotential ( )
virtual

Reimplemented in SmallStrainJ2Plasticity, and SmallStrainParaboloidalPlasticity.

Examples
SmallStrainPlasticity.hpp.

Definition at line 591 of file SmallStrainPlasticity.hpp.

591 {
592 PetscFunctionBegin;
593 SETERRQ(
594 PETSC_COMM_SELF,
596 "Not implemented (overload class and implement this function)"
597 );
598 PetscFunctionReturn(0);
599 }
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32

◆ freeHemholtzEnergy()

virtual PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::freeHemholtzEnergy ( )
virtual

Reimplemented in SmallStrainJ2Plasticity, and SmallStrainParaboloidalPlasticity.

Examples
SmallStrainPlasticity.hpp.

Definition at line 571 of file SmallStrainPlasticity.hpp.

571 {
572 PetscFunctionBegin;
573 SETERRQ(
574 PETSC_COMM_SELF,
576 "Not implemented (overload class and implement this function)"
577 );
578 PetscFunctionReturn(0);
579 }

◆ pLayH()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::pLayH ( )
virtual
Examples
SmallStrainPlasticity.hpp.

Definition at line 855 of file SmallStrainPlasticity.cpp.

855 {
856 PetscFunctionBegin;
857 int r;
858 int adloc_return_value = 0;
859 r = ::function(tAgs[2],1,activeVariablesYH.size(),&activeVariablesYH[0],&h);
860 if(r<adloc_return_value) {
861 SETERRQ(PETSC_COMM_SELF,MOFEM_OPERATION_UNSUCCESSFUL,"ADOL-C function evaluation with error");
862 }
863 gradientH.resize(activeVariablesYH.size());
864 r = gradient(tAgs[2],activeVariablesYH.size(),&activeVariablesYH[0],&gradientH[0]);
865 if(r<adloc_return_value) {
866 SETERRQ(PETSC_COMM_SELF,MOFEM_OPERATION_UNSUCCESSFUL,"ADOL-C function evaluation with error");
867 }
868 hessianH.resize(6+internalVariables.size(),6+internalVariables.size(),false);
869 hessianH.clear();
870 vector<double*> hessian_h(activeVariablesYH.size());
871 for(int dd = 0;dd<activeVariablesYH.size();dd++) {
872 hessian_h[dd] = &hessianH(dd,0);
873 }
874 r = hessian(tAgs[2],activeVariablesYH.size(),&activeVariablesYH[0],&hessian_h[0]);
875 if(r<adloc_return_value) {
876 SETERRQ(PETSC_COMM_SELF,MOFEM_OPERATION_UNSUCCESSFUL,"ADOL-C function evaluation with error");
877 }
878 partialHSigma = VectorAdaptor(6,ublas::shallow_array_adaptor<double>(6,&gradientH[0]));
880 internalVariables.size(),
881 ublas::shallow_array_adaptor<double>(
882 internalVariables.size(),
883 &gradientH[6]
884 )
885 );
886 partialHFlux *= -1;
887 partial2HSigma.resize(6,6,false);
888 for(int ii = 0;ii<6;ii++) {
889 for(int jj = 0;jj<=ii;jj++) {
890 partial2HSigma(ii,jj) = hessianH(ii,jj);
891 }
892 }
893 partial2HFlux.resize(internalVariables.size(),internalVariables.size(),false);
894 for(unsigned int ii = 0;ii<internalVariables.size();ii++) {
895 for(unsigned int jj = 0;jj<=ii;jj++) {
896 partial2HFlux(ii,jj) = -hessianH(6+ii,6+jj);
897 }
898 }
899 partial2HSigmaFlux.resize(6,internalVariables.size(),false);
900 partial2HSigmaFlux.clear();
901 for(int ii = 0;ii<6;ii++) {
902 for(unsigned int jj = 0;jj<internalVariables.size();jj++) {
903 partial2HSigmaFlux(ii,jj) = -hessianH(6+jj,ii);
904 }
905 }
906 PetscFunctionReturn(0);
907}
@ MOFEM_OPERATION_UNSUCCESSFUL
Definition: definitions.h:34
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition: Types.hpp:115
const double r
rate factor

◆ pLayW()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::pLayW ( )
virtual
Examples
SmallStrainPlasticity.hpp.

Definition at line 747 of file SmallStrainPlasticity.cpp.

747 {
748 PetscFunctionBegin;
749 int r;
750 int adloc_return_value = 0;
751 r = ::function(tAgs[0],1,activeVariablesW.size(),&activeVariablesW[0],&w);
752 if(r<adloc_return_value) {
753 SETERRQ(PETSC_COMM_SELF,MOFEM_OPERATION_UNSUCCESSFUL,"ADOL-C function evaluation with error");
754 }
755 gradientW.resize(activeVariablesW.size(),false);
756 r = gradient(tAgs[0],activeVariablesW.size(),&activeVariablesW[0],&gradientW[0]);
757 if(r<adloc_return_value) {
758 SETERRQ(PETSC_COMM_SELF,MOFEM_OPERATION_UNSUCCESSFUL,"ADOL-C function evaluation with error");
759 }
760 hessianW.resize(activeVariablesW.size(),activeVariablesW.size(),false);
761 hessianW.clear();
762 vector<double*> hessian_w(activeVariablesW.size());
763 for(unsigned int dd = 0;dd<activeVariablesW.size();dd++) {
764 hessian_w[dd] = &hessianW(dd,0);
765 }
766 r = hessian(tAgs[0],activeVariablesW.size(),&activeVariablesW[0],&hessian_w[0]);
767 if(r<adloc_return_value) {
768 SETERRQ(PETSC_COMM_SELF,MOFEM_OPERATION_UNSUCCESSFUL,"ADOL-C function evaluation with error");
769 }
770 sTress = VectorAdaptor(6,ublas::shallow_array_adaptor<double>(6,&gradientW[0]));
772 internalVariables.size(),
773 ublas::shallow_array_adaptor<double>(
774 internalVariables.size(),
775 &gradientW[12]
776 )
777 );
778 C.resize(6,6,false);
779 for(int ii = 0;ii<6;ii++) {
780 for(int jj = 0;jj<=ii;jj++) {
781 C(ii,jj) = hessianW(ii,jj);
782 }
783 }
784 partialWStrainPlasticStrain.resize(6,6,false);
785 for(int ii = 0;ii<6;ii++) {
786 for(int jj = 0;jj<6;jj++) {
787 partialWStrainPlasticStrain(ii,jj) = hessianW(6+jj,ii);
788 }
789 }
790 D.resize(internalVariables.size(),internalVariables.size(),false);
791 for(unsigned int ii = 0;ii<internalVariables.size();ii++) {
792 for(unsigned int jj = 0;jj<=ii;jj++) {
793 D(ii,jj) = hessianW(12+ii,12+jj);
794 }
795 }
796 PetscFunctionReturn(0);
797}

◆ pLayW_NoHessian()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::pLayW_NoHessian ( )
virtual
Examples
SmallStrainPlasticity.hpp.

Definition at line 798 of file SmallStrainPlasticity.cpp.

798 {
799 PetscFunctionBegin;
800 int r;
801 int adloc_return_value = 0;
802 r = ::function(tAgs[0],1,activeVariablesW.size(),&activeVariablesW[0],&w);
803 if(r<adloc_return_value) {
804 SETERRQ(PETSC_COMM_SELF,MOFEM_OPERATION_UNSUCCESSFUL,"ADOL-C function evaluation with error");
805 }
806 gradientW.resize(activeVariablesW.size(),false);
807 r = gradient(tAgs[0],activeVariablesW.size(),&activeVariablesW[0],&gradientW[0]);
808 if(r<adloc_return_value) {
809 SETERRQ(PETSC_COMM_SELF,MOFEM_OPERATION_UNSUCCESSFUL,"ADOL-C function evaluation with error");
810 }
811 sTress = VectorAdaptor(6,ublas::shallow_array_adaptor<double>(6,&gradientW[0]));
813 internalVariables.size(),
814 ublas::shallow_array_adaptor<double>(
815 internalVariables.size(),
816 &gradientW[12]
817 )
818 );
819 PetscFunctionReturn(0);
820}

◆ pLayY()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::pLayY ( )
virtual
Examples
SmallStrainPlasticity.hpp.

Definition at line 822 of file SmallStrainPlasticity.cpp.

822 {
823 PetscFunctionBegin;
824 int r;
825 int adloc_return_value = 0;
826 r = ::function(tAgs[1],1,activeVariablesYH.size(),&activeVariablesYH[0],&y);
827 if(r<adloc_return_value) {
828 SETERRQ(PETSC_COMM_SELF,MOFEM_OPERATION_UNSUCCESSFUL,"ADOL-C function evaluation with error");
829 }
830 gradientY.resize(activeVariablesYH.size());
831 r = gradient(tAgs[1],activeVariablesYH.size(),&activeVariablesYH[0],&gradientY[0]);
832 if(r<adloc_return_value) {
833 SETERRQ(PETSC_COMM_SELF,MOFEM_OPERATION_UNSUCCESSFUL,"ADOL-C function evaluation with error");
834 }
835 partialYSigma = VectorAdaptor(6,ublas::shallow_array_adaptor<double>(6,&gradientY[0]));
837 internalVariables.size(),
838 ublas::shallow_array_adaptor<double>(
839 internalVariables.size(),
840 &gradientY[6]
841 )
842 );
843 PetscFunctionReturn(0);
844}

◆ pLayY_NoGradient()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::pLayY_NoGradient ( )
virtual
Examples
SmallStrainPlasticity.hpp.

Definition at line 845 of file SmallStrainPlasticity.cpp.

845 {
846 PetscFunctionBegin;
847 int r;
848 int adloc_return_value = 0;
849 r = ::function(tAgs[1],1,activeVariablesYH.size(),&activeVariablesYH[0],&y);
850 if(r<adloc_return_value) {
851 SETERRQ(PETSC_COMM_SELF,MOFEM_OPERATION_UNSUCCESSFUL,"ADOL-C function evaluation with error");
852 }
853 PetscFunctionReturn(0);
854}

◆ rEcordH()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::rEcordH ( )
virtual
Examples
SmallStrainPlasticity.hpp.

Definition at line 725 of file SmallStrainPlasticity.cpp.

725 {
726 PetscFunctionBegin;
727 trace_on(tAgs[2]);
728 {
729 //active variables
730 a_sTress.resize(6,false);
731 for(int ii = 0;ii<6;ii++) {
732 a_sTress[ii] <<= sTress[ii];
733 }
734 a_internalFluxes.resize(internalVariables.size(),false);
735 for(unsigned int ii = 0;ii<internalVariables.size();ii++) {
737 }
738 //evaluete functions
739 ierr = flowPotential(); CHKERRQ(ierr);
740 //return variables
741 a_h >>= h;
742 }
743 trace_off();
744 PetscFunctionReturn(0);
745}

◆ rEcordW()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::rEcordW ( )
virtual
Examples
SmallStrainPlasticity.hpp.

Definition at line 677 of file SmallStrainPlasticity.cpp.

677 {
678 PetscFunctionBegin;
679 trace_on(tAgs[0]);
680 {
681 //active variables
682 a_sTrain.resize(6,false);
683 for(int ii = 0;ii<6;ii++) {
684 a_sTrain[ii] <<= sTrain[ii];
685 }
686 a_plasticStrain.resize(6,false);
687 for(int ii = 0;ii<6;ii++) {
688 a_plasticStrain[ii] <<= plasticStrain[ii];
689 }
690 a_internalVariables.resize(internalVariables.size(),false);
691 for(unsigned int ii = 0;ii<internalVariables.size();ii++) {
693 }
694 //evaluete functions
695 ierr = freeHemholtzEnergy(); CHKERRQ(ierr);
696 //return variables
697 a_w >>= w;
698 }
699 trace_off();
700 PetscFunctionReturn(0);
701}

◆ rEcordY()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::rEcordY ( )
virtual
Examples
SmallStrainPlasticity.hpp.

Definition at line 703 of file SmallStrainPlasticity.cpp.

703 {
704 PetscFunctionBegin;
705 trace_on(tAgs[1]);
706 {
707 //active variables
708 a_sTress.resize(6,false);
709 for(int ii = 0;ii<6;ii++) {
710 a_sTress[ii] <<= sTress[ii];
711 }
712 a_internalFluxes.resize(internalVariables.size(),false);
713 for(unsigned int ii = 0;ii<internalVariables.size();ii++) {
715 }
716 //evaluete functions
717 ierr = yieldFunction(); CHKERRQ(ierr);
718 //return variables
719 a_y >>= y;
720 }
721 trace_off();
722 PetscFunctionReturn(0);
723}

◆ setActiveVariablesW()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::setActiveVariablesW ( )
virtual
Examples
SmallStrainPlasticity.hpp.

Definition at line 650 of file SmallStrainPlasticity.cpp.

650 {
651 PetscFunctionBegin;
652 activeVariablesW.resize(12+internalVariables.size(),false);
653 for(int ii = 0;ii<6;ii++) {
654 activeVariablesW[ii] = sTrain[ii];
655 }
656 for(int ii = 0;ii<6;ii++) {
658 }
659 for(unsigned int ii = 0;ii<internalVariables.size();ii++) {
661 }
662 PetscFunctionReturn(0);
663}

◆ setActiveVariablesYH()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::setActiveVariablesYH ( )
virtual
Examples
SmallStrainPlasticity.hpp.

Definition at line 665 of file SmallStrainPlasticity.cpp.

665 {
666 PetscFunctionBegin;
667 activeVariablesYH.resize(6+internalVariables.size(),false);
668 for(int ii = 0;ii<6;ii++) {
669 activeVariablesYH[ii] = sTress[ii];
670 }
671 for(unsigned int ii = 0;ii<internalVariables.size();ii++) {
673 }
674 PetscFunctionReturn(0);
675}

◆ snesCreate()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::snesCreate ( )
Examples
SmallStrainPlasticity.hpp.

Definition at line 1012 of file SmallStrainPlasticity.cpp.

1012 {
1013 PetscFunctionBegin;
1014 ierr = SNESCreate(PETSC_COMM_SELF,&sNes); CHKERRQ(ierr);
1015 ierr = SNESSetFunction(sNes,R,SmallStrainPlasticityfRes,(void*)this); CHKERRQ(ierr);
1016 ierr = SNESSetJacobian(sNes,A,A,SmallStrainPlasticityfJac,(void*)this); CHKERRQ(ierr);
1017 ierr = SNESGetKSP(sNes,&kSp); CHKERRQ(ierr);
1018 ierr = KSPGetPC(kSp,&pC); CHKERRQ(ierr);
1019 ierr = KSPSetType(kSp,KSPPREONLY); CHKERRQ(ierr);
1020 ierr = PCSetType(pC,PCLU); CHKERRQ(ierr);
1021 ierr = SNESSetTolerances(sNes,tOl,1e-12,0,20,1000); CHKERRQ(ierr);
1022 ierr = SNESGetLineSearch(sNes,&lineSearch); CHKERRQ(ierr); CHKERRQ(ierr);
1023 //ierr = SNESLineSearchSetType(lineSearch,SNESLINESEARCHBASIC); CHKERRQ(ierr);
1024 ierr = SNESLineSearchSetType(lineSearch,SNESLINESEARCHBT); CHKERRQ(ierr);
1025
1026 PetscFunctionReturn(0);
1027}
friend PetscErrorCode SmallStrainPlasticityfJac(SNES, Vec, Mat, Mat, void *ctx)
friend PetscErrorCode SmallStrainPlasticityfRes(SNES, Vec, Vec, void *ctx)

◆ snesDestroy()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::snesDestroy ( )
Examples
SmallStrainPlasticity.hpp.

Definition at line 1029 of file SmallStrainPlasticity.cpp.

1029 {
1030 PetscFunctionBegin;
1031 ierr = SNESDestroy(&sNes); CHKERRQ(ierr);
1032 PetscFunctionReturn(0);
1033}

◆ solveColasetProjection()

PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::solveColasetProjection ( )
Examples
SmallStrainPlasticity.hpp.

Definition at line 1035 of file SmallStrainPlasticity.cpp.

1035 {
1036 PetscFunctionBegin;
1037 {
1038 double *array;
1039 ierr = VecGetArray(Chi,&array); CHKERRQ(ierr);
1040 for(int ii = 0;ii<6;ii++) {
1041 array[ii] = plasticStrain[ii];
1042 }
1043 int nb_internal_variables = internalVariables.size();
1044 for(unsigned int ii = 0;ii<nb_internal_variables;ii++) {
1045 array[6+ii] = internalVariables[ii];
1046 }
1047 deltaGamma = 0;
1048 ierr = VecRestoreArray(Chi,&array); CHKERRQ(ierr);
1049 }
1050 ierr = SNESSolve(sNes,PETSC_NULL,Chi); CHKERRQ(ierr);
1051 SNESConvergedReason reason;
1052 ierr = SNESGetConvergedReason(sNes,&reason); CHKERRQ(ierr);
1053 if(reason < 0) {
1054 int its;
1055 ierr = SNESGetIterationNumber(sNes,&its); CHKERRQ(ierr);
1056// ierr = PetscPrintf(
1057// PETSC_COMM_SELF,
1058// "** WARNING Plasticity Closet Point Projection - number of Newton iterations = %D\n",
1059// its
1060// ); CHKERRQ(ierr);
1061 noalias(plasticStrain) = plasticStrain0;
1063 deltaGamma = 0;
1064 } else {
1065 double *array;
1066 ierr = VecGetArray(Chi,&array); CHKERRQ(ierr);
1067 for(int ii = 0;ii<6;ii++) {
1068 plasticStrain[ii] = array[ii];
1069 }
1070 int nb_internal_variables = internalVariables.size();
1071 for(unsigned int ii = 0;ii<nb_internal_variables;ii++) {
1072 internalVariables[ii] = array[6+ii];
1073 }
1074 deltaGamma = array[6+nb_internal_variables];
1075 ierr = VecRestoreArray(Chi,&array); CHKERRQ(ierr);
1076 }
1077 PetscFunctionReturn(0);
1078}

◆ yieldFunction()

virtual PetscErrorCode SmallStrainPlasticity::ClosestPointProjection::yieldFunction ( )
virtual

Reimplemented in SmallStrainJ2Plasticity, and SmallStrainParaboloidalPlasticity.

Examples
SmallStrainPlasticity.hpp.

Definition at line 581 of file SmallStrainPlasticity.hpp.

581 {
582 PetscFunctionBegin;
583 SETERRQ(
584 PETSC_COMM_SELF,
586 "Not implemented (overload class and implement this function)"
587 );
588 PetscFunctionReturn(0);
589 }

Friends And Related Function Documentation

◆ SmallStrainPlasticityfJac

PetscErrorCode SmallStrainPlasticityfJac ( SNES  ,
Vec  ,
Mat  ,
Mat  ,
void *  ctx 
)
friend
Examples
SmallStrainPlasticity.hpp.

Definition at line 1192 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}

◆ SmallStrainPlasticityfRes

PetscErrorCode SmallStrainPlasticityfRes ( SNES  ,
Vec  ,
Vec  ,
void *  ctx 
)
friend
Examples
SmallStrainPlasticity.hpp.

Definition at line 1139 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}

Member Data Documentation

◆ A

Mat SmallStrainPlasticity::ClosestPointProjection::A
Examples
SmallStrainPlasticity.hpp.

Definition at line 540 of file SmallStrainPlasticity.hpp.

◆ a_h

adouble SmallStrainPlasticity::ClosestPointProjection::a_h

◆ a_internalFluxes

ublas::vector<adouble> SmallStrainPlasticity::ClosestPointProjection::a_internalFluxes

◆ a_internalVariables

ublas::vector<adouble> SmallStrainPlasticity::ClosestPointProjection::a_internalVariables

◆ a_plasticStrain

ublas::vector<adouble> SmallStrainPlasticity::ClosestPointProjection::a_plasticStrain

◆ a_sTrain

ublas::vector<adouble> SmallStrainPlasticity::ClosestPointProjection::a_sTrain

◆ a_sTress

ublas::vector<adouble> SmallStrainPlasticity::ClosestPointProjection::a_sTress

◆ a_w

adouble SmallStrainPlasticity::ClosestPointProjection::a_w

◆ a_y

adouble SmallStrainPlasticity::ClosestPointProjection::a_y

◆ activeVariablesW

VectorDouble SmallStrainPlasticity::ClosestPointProjection::activeVariablesW
Examples
SmallStrainPlasticity.hpp.

Definition at line 510 of file SmallStrainPlasticity.hpp.

◆ activeVariablesYH

VectorDouble SmallStrainPlasticity::ClosestPointProjection::activeVariablesYH
Examples
SmallStrainPlasticity.hpp.

Definition at line 513 of file SmallStrainPlasticity.hpp.

◆ C

ublas::symmetric_matrix<double,ublas::lower> SmallStrainPlasticity::ClosestPointProjection::C
Examples
SmallStrainPlasticity.hpp.

Definition at line 496 of file SmallStrainPlasticity.hpp.

◆ Cep

MatrixDouble SmallStrainPlasticity::ClosestPointProjection::Cep
Examples
SmallStrainPlasticity.hpp.

Definition at line 566 of file SmallStrainPlasticity.hpp.

◆ Chi

Vec SmallStrainPlasticity::ClosestPointProjection::Chi
Examples
SmallStrainPlasticity.hpp.

Definition at line 542 of file SmallStrainPlasticity.hpp.

◆ Cp

MatrixDouble SmallStrainPlasticity::ClosestPointProjection::Cp
Examples
SmallStrainPlasticity.hpp.

Definition at line 566 of file SmallStrainPlasticity.hpp.

◆ D

ublas::symmetric_matrix<double,ublas::lower> SmallStrainPlasticity::ClosestPointProjection::D
Examples
SmallStrainPlasticity.hpp.

Definition at line 496 of file SmallStrainPlasticity.hpp.

◆ dataA

ublas::matrix<double,ublas::column_major> SmallStrainPlasticity::ClosestPointProjection::dataA
Examples
SmallStrainPlasticity.hpp.

Definition at line 541 of file SmallStrainPlasticity.hpp.

◆ dChi

Vec SmallStrainPlasticity::ClosestPointProjection::dChi
Examples
SmallStrainPlasticity.hpp.

Definition at line 542 of file SmallStrainPlasticity.hpp.

◆ deltaGamma

double SmallStrainPlasticity::ClosestPointProjection::deltaGamma
Examples
SmallStrainPlasticity.hpp.

Definition at line 490 of file SmallStrainPlasticity.hpp.

◆ Ep

MatrixDouble SmallStrainPlasticity::ClosestPointProjection::Ep
Examples
SmallStrainPlasticity.hpp.

Definition at line 566 of file SmallStrainPlasticity.hpp.

◆ gG

int SmallStrainPlasticity::ClosestPointProjection::gG
Examples
SmallStrainPlasticity.hpp.

Definition at line 492 of file SmallStrainPlasticity.hpp.

◆ gradientH

VectorDouble SmallStrainPlasticity::ClosestPointProjection::gradientH
Examples
SmallStrainPlasticity.hpp.

Definition at line 536 of file SmallStrainPlasticity.hpp.

◆ gradientW

VectorDouble SmallStrainPlasticity::ClosestPointProjection::gradientW
Examples
SmallStrainPlasticity.hpp.

Definition at line 527 of file SmallStrainPlasticity.hpp.

◆ gradientY

VectorDouble SmallStrainPlasticity::ClosestPointProjection::gradientY
Examples
SmallStrainPlasticity.hpp.

Definition at line 532 of file SmallStrainPlasticity.hpp.

◆ h

double SmallStrainPlasticity::ClosestPointProjection::h
Examples
SmallStrainPlasticity.hpp.

Definition at line 516 of file SmallStrainPlasticity.hpp.

◆ hessianH

MatrixDouble SmallStrainPlasticity::ClosestPointProjection::hessianH
Examples
SmallStrainPlasticity.hpp.

Definition at line 537 of file SmallStrainPlasticity.hpp.

◆ hessianW

MatrixDouble SmallStrainPlasticity::ClosestPointProjection::hessianW
Examples
SmallStrainPlasticity.hpp.

Definition at line 528 of file SmallStrainPlasticity.hpp.

◆ internalFluxes

VectorAdaptor SmallStrainPlasticity::ClosestPointProjection::internalFluxes
Examples
SmallStrainPlasticity.hpp.

Definition at line 495 of file SmallStrainPlasticity.hpp.

◆ internalVariables

VectorDouble SmallStrainPlasticity::ClosestPointProjection::internalVariables

◆ internalVariables0

VectorDouble SmallStrainPlasticity::ClosestPointProjection::internalVariables0
Examples
SmallStrainPlasticity.hpp.

Definition at line 487 of file SmallStrainPlasticity.hpp.

◆ kSp

KSP SmallStrainPlasticity::ClosestPointProjection::kSp
Examples
SmallStrainPlasticity.hpp.

Definition at line 557 of file SmallStrainPlasticity.hpp.

◆ lineSearch

SNESLineSearch SmallStrainPlasticity::ClosestPointProjection::lineSearch
Examples
SmallStrainPlasticity.hpp.

Definition at line 559 of file SmallStrainPlasticity.hpp.

◆ Lp

MatrixDouble SmallStrainPlasticity::ClosestPointProjection::Lp
Examples
SmallStrainPlasticity.hpp.

Definition at line 566 of file SmallStrainPlasticity.hpp.

◆ partial2HFlux

ublas::symmetric_matrix<double,ublas::lower> SmallStrainPlasticity::ClosestPointProjection::partial2HFlux
Examples
SmallStrainPlasticity.hpp.

Definition at line 503 of file SmallStrainPlasticity.hpp.

◆ partial2HSigma

ublas::symmetric_matrix<double,ublas::lower> SmallStrainPlasticity::ClosestPointProjection::partial2HSigma
Examples
SmallStrainPlasticity.hpp.

Definition at line 502 of file SmallStrainPlasticity.hpp.

◆ partial2HSigmaFlux

MatrixDouble SmallStrainPlasticity::ClosestPointProjection::partial2HSigmaFlux
Examples
SmallStrainPlasticity.hpp.

Definition at line 504 of file SmallStrainPlasticity.hpp.

◆ partialHFlux

VectorAdaptor SmallStrainPlasticity::ClosestPointProjection::partialHFlux
Examples
SmallStrainPlasticity.hpp.

Definition at line 501 of file SmallStrainPlasticity.hpp.

◆ partialHSigma

VectorAdaptor SmallStrainPlasticity::ClosestPointProjection::partialHSigma
Examples
SmallStrainPlasticity.hpp.

Definition at line 500 of file SmallStrainPlasticity.hpp.

◆ partialWStrainPlasticStrain

MatrixDouble SmallStrainPlasticity::ClosestPointProjection::partialWStrainPlasticStrain
Examples
SmallStrainPlasticity.hpp.

Definition at line 497 of file SmallStrainPlasticity.hpp.

◆ partialYFlux

VectorAdaptor SmallStrainPlasticity::ClosestPointProjection::partialYFlux
Examples
SmallStrainPlasticity.hpp.

Definition at line 499 of file SmallStrainPlasticity.hpp.

◆ partialYSigma

VectorAdaptor SmallStrainPlasticity::ClosestPointProjection::partialYSigma
Examples
SmallStrainPlasticity.hpp.

Definition at line 498 of file SmallStrainPlasticity.hpp.

◆ pC

PC SmallStrainPlasticity::ClosestPointProjection::pC
Examples
SmallStrainPlasticity.hpp.

Definition at line 558 of file SmallStrainPlasticity.hpp.

◆ plasticStrain

VectorDouble SmallStrainPlasticity::ClosestPointProjection::plasticStrain
Examples
SmallStrainPlasticity.hpp.

Definition at line 489 of file SmallStrainPlasticity.hpp.

◆ plasticStrain0

VectorDouble SmallStrainPlasticity::ClosestPointProjection::plasticStrain0
Examples
SmallStrainPlasticity.hpp.

Definition at line 488 of file SmallStrainPlasticity.hpp.

◆ R

Vec SmallStrainPlasticity::ClosestPointProjection::R
Examples
SmallStrainPlasticity.hpp.

Definition at line 542 of file SmallStrainPlasticity.hpp.

◆ sNes

SNES SmallStrainPlasticity::ClosestPointProjection::sNes
Examples
SmallStrainPlasticity.hpp.

Definition at line 556 of file SmallStrainPlasticity.hpp.

◆ sTrain

VectorDouble SmallStrainPlasticity::ClosestPointProjection::sTrain
Examples
SmallStrainPlasticity.hpp.

Definition at line 485 of file SmallStrainPlasticity.hpp.

◆ sTress

VectorAdaptor SmallStrainPlasticity::ClosestPointProjection::sTress
Examples
SmallStrainPlasticity.hpp.

Definition at line 494 of file SmallStrainPlasticity.hpp.

◆ tAgs

vector<int> SmallStrainPlasticity::ClosestPointProjection::tAgs
Examples
SmallStrainPlasticity.hpp.

Definition at line 506 of file SmallStrainPlasticity.hpp.

◆ tOl

double SmallStrainPlasticity::ClosestPointProjection::tOl
Examples
SmallStrainPlasticity.hpp.

Definition at line 491 of file SmallStrainPlasticity.hpp.

◆ w

double SmallStrainPlasticity::ClosestPointProjection::w
Examples
SmallStrainPlasticity.hpp.

Definition at line 516 of file SmallStrainPlasticity.hpp.

◆ y

double SmallStrainPlasticity::ClosestPointProjection::y
Examples
SmallStrainPlasticity.hpp.

Definition at line 516 of file SmallStrainPlasticity.hpp.

◆ Yp

VectorDouble SmallStrainPlasticity::ClosestPointProjection::Yp
Examples
SmallStrainPlasticity.hpp.

Definition at line 567 of file SmallStrainPlasticity.hpp.


The documentation for this struct was generated from the following files: