 |
Ice Sheet System Model
4.18
Code documentation
|
#include <MasstransportAnalysis.h>
|
void | CreateConstraints (Constraints *constraints, IoModel *iomodel) |
|
void | CreateLoads (Loads *loads, IoModel *iomodel) |
|
void | CreateNodes (Nodes *nodes, IoModel *iomodel, bool isamr=false) |
|
int | DofsPerNode (int **doflist, int domaintype, int approximation) |
|
void | UpdateElements (Elements *elements, Inputs2 *inputs2, IoModel *iomodel, int analysis_counter, int analysis_type) |
|
void | UpdateParameters (Parameters *parameters, IoModel *iomodel, int solution_enum, int analysis_enum) |
|
void | Core (FemModel *femmodel) |
|
ElementVector * | CreateDVector (Element *element) |
|
ElementMatrix * | CreateJacobianMatrix (Element *element) |
|
ElementMatrix * | CreateKMatrix (Element *element) |
|
ElementMatrix * | CreateKMatrixCG (Element *element) |
|
ElementMatrix * | CreateKMatrixDG (Element *element) |
|
ElementVector * | CreatePVector (Element *element) |
|
ElementVector * | CreatePVectorCG (Element *element) |
|
ElementVector * | CreatePVectorDG (Element *element) |
|
void | GetSolutionFromInputs (Vector< IssmDouble > *solution, Element *element) |
|
void | GradientJ (Vector< IssmDouble > *gradient, Element *element, int control_type, int control_index) |
|
void | InputUpdateFromSolution (IssmDouble *solution, Element *element) |
|
void | UpdateConstraints (FemModel *femmodel) |
|
ElementMatrix * | CreateFctKMatrix (Element *element) |
|
ElementMatrix * | CreateMassMatrix (Element *element) |
|
ElementVector * | CreateFctPVector (Element *element) |
|
void | FctKMatrix (Matrix< IssmDouble > **pKff, Matrix< IssmDouble > **pKfs, FemModel *femmodel) |
|
void | FctPVector (Vector< IssmDouble > **ppf, FemModel *femmodel) |
|
void | MassMatrix (Matrix< IssmDouble > **pMff, FemModel *femmodel) |
|
virtual | ~Analysis () |
|
Definition at line 11 of file MasstransportAnalysis.h.
◆ CreateConstraints()
void MasstransportAnalysis::CreateConstraints |
( |
Constraints * |
constraints, |
|
|
IoModel * |
iomodel |
|
) |
| |
|
virtual |
◆ CreateLoads()
void MasstransportAnalysis::CreateLoads |
( |
Loads * |
loads, |
|
|
IoModel * |
iomodel |
|
) |
| |
|
virtual |
Implements Analysis.
Definition at line 26 of file MasstransportAnalysis.cpp.
32 int numvertex_pairing;
35 iomodel->
FindConstant(&stabilization,
"md.masstransport.stabilization");
42 iomodel->
FetchData(1,
"md.geometry.thickness");
48 int element=iomodel->
faces[4*i+2]-1;
56 iomodel->
DeleteData(1,
"md.geometry.thickness");
62 iomodel->
FetchData(&vertex_pairing,&numvertex_pairing,NULL,
"md.masstransport.vertex_pairing");
65 for(
int i=0;i<numvertex_pairing;i++){
67 if(iomodel->
my_vertices[reCast<int>(vertex_pairing[2*i+0])-1]){
74 if(!(reCast<bool>(nodeonbase[reCast<int>(vertex_pairing[2*i+0])-1])) || !(reCast<bool>(nodeonbase[reCast<int>(vertex_pairing[2*i+1])-1])))
continue;
78 penpair_ids[0]=reCast<int>(vertex_pairing[2*i+0]);
79 penpair_ids[1]=reCast<int>(vertex_pairing[2*i+1]);
88 iomodel->
DeleteData(vertex_pairing,
"md.masstransport.vertex_pairing");
89 iomodel->
DeleteData(nodeonbase,
"md.mesh.vertexonbase");
◆ CreateNodes()
void MasstransportAnalysis::CreateNodes |
( |
Nodes * |
nodes, |
|
|
IoModel * |
iomodel, |
|
|
bool |
isamr = false |
|
) |
| |
|
virtual |
◆ DofsPerNode()
int MasstransportAnalysis::DofsPerNode |
( |
int ** |
doflist, |
|
|
int |
domaintype, |
|
|
int |
approximation |
|
) |
| |
|
virtual |
◆ UpdateElements()
void MasstransportAnalysis::UpdateElements |
( |
Elements * |
elements, |
|
|
Inputs2 * |
inputs2, |
|
|
IoModel * |
iomodel, |
|
|
int |
analysis_counter, |
|
|
int |
analysis_type |
|
) |
| |
|
virtual |
Implements Analysis.
Definition at line 113 of file MasstransportAnalysis.cpp.
115 int stabilization,finiteelement;
116 bool dakota_analysis;
117 bool isgroundingline;
119 bool isoceancoupling;
123 iomodel->
FindConstant(&stabilization,
"md.masstransport.stabilization");
124 iomodel->
FindConstant(&dakota_analysis,
"md.qmu.isdakota");
125 iomodel->
FindConstant(&isgroundingline,
"md.transient.isgroundingline");
126 iomodel->
FindConstant(&ismovingfront,
"md.transient.ismovingfront");
127 iomodel->
FindConstant(&isoceancoupling,
"md.transient.isoceancoupling");
132 if(stabilization==3){
142 element->
Update(inputs2,i,iomodel,analysis_counter,analysis_type,finiteelement);
163 int basalforcing_model;
164 iomodel->
FindConstant(&basalforcing_model,
"md.basalforcings.model");
165 switch(basalforcing_model){
188 IssmDouble** array3d = NULL;
int* Ms = NULL;
int* Ns = NULL;
int K;
189 iomodel->
FetchData(&array3d,&Ms,&Ns,&K,
"md.basalforcings.tf");
190 if(!array3d)
_error_(
"md.basalforcings.tf not found in binary file");
191 for(
int i=0;i<elements->
Size();i++){
194 for(
int kk=0;kk<K;kk++){
198 xDelete<int>(Ms); xDelete<int>(Ns);
199 for(
int i=0;i<K;i++) xDelete<IssmDouble>(array3d[i]);
200 xDelete<IssmDouble*>(array3d);
214 if(stabilization==3){
217 if(stabilization==4){
◆ UpdateParameters()
void MasstransportAnalysis::UpdateParameters |
( |
Parameters * |
parameters, |
|
|
IoModel * |
iomodel, |
|
|
int |
solution_enum, |
|
|
int |
analysis_enum |
|
) |
| |
|
virtual |
Implements Analysis.
Definition at line 227 of file MasstransportAnalysis.cpp.
230 char** requestedoutputs = NULL;
239 iomodel->
FindConstant(&requestedoutputs,&numoutputs,
"md.masstransport.requested_outputs");
242 iomodel->
DeleteData(&requestedoutputs,numoutputs,
"md.masstransport.requested_outputs");
◆ Core()
void MasstransportAnalysis::Core |
( |
FemModel * |
femmodel | ) |
|
|
virtual |
◆ CreateDVector()
◆ CreateJacobianMatrix()
◆ CreateKMatrix()
◆ CreateKMatrixCG()
Definition at line 283 of file MasstransportAnalysis.cpp.
312 IssmDouble* basis = xNew<IssmDouble>(numnodes);
313 IssmDouble* dbasis = xNew<IssmDouble>(dim*numnodes);
321 Input2* vyaverage_input=NULL;
330 for(
int ig=gauss->
begin();ig<gauss->end();ig++){
338 D_scalar=gauss->
weight*Jdet;
339 for(
int i=0;i<numnodes;i++){
340 for(
int j=0;j<numnodes;j++){
341 Ke->
values[i*numnodes+j] += D_scalar*basis[i]*basis[j];
348 D_scalar=dt*gauss->
weight*Jdet;
354 for(
int i=0;i<numnodes;i++){
355 for(
int j=0;j<numnodes;j++){
357 Ke->
values[i*numnodes+j] += D_scalar*basis[i]*basis[j]*(dvxdx+dvydy);
359 Ke->
values[i*numnodes+j] += D_scalar*basis[i]*(vx*dbasis[0*numnodes+j] + vy*dbasis[1*numnodes+j]);
365 for(
int i=0;i<numnodes;i++){
366 for(
int j=0;j<numnodes;j++){
367 Ke->
values[i*numnodes+j] += D_scalar*dvxdx*dbasis[0*numnodes+i]*dbasis[0*numnodes+j];
368 Ke->
values[i*numnodes+j] += D_scalar*vx*dbasis[0*numnodes+j]*basis[i];
373 for(
int i=0;i<4;i++) D[i]=0.;
374 switch(stabilization){
382 D[0*dim+0]=h/2.0*fabs(vx);
383 if(dim==2) D[1*dim+1]=h/2.0*fabs(vy);
393 vel=sqrt(vx*vx+vy*vy)+1.e-8;
399 if(dim!=2)
_error_(
"Stabilization "<<stabilization<<
" not supported yet for dim != 2");
402 vel=sqrt(vx*vx+vy*vy)+1.e-8;
408 _error_(
"Stabilization "<<stabilization<<
" not supported yet");
410 if(stabilization==1){
412 if(dim==1) D[0]=D_scalar*D[0];
414 D[0*dim+0]=D_scalar*D[0*dim+0];
415 D[1*dim+0]=D_scalar*D[1*dim+0];
416 D[0*dim+1]=D_scalar*D[0*dim+1];
417 D[1*dim+1]=D_scalar*D[1*dim+1];
421 for(
int i=0;i<numnodes;i++){
422 for(
int j=0;j<numnodes;j++){
423 Ke->
values[i*numnodes+j] += (
424 dbasis[0*numnodes+i] *(D[0*dim+0]*dbasis[0*numnodes+j] + D[0*dim+1]*dbasis[1*numnodes+j]) +
425 dbasis[1*numnodes+i] *(D[1*dim+0]*dbasis[0*numnodes+j] + D[1*dim+1]*dbasis[1*numnodes+j])
431 for(
int i=0;i<numnodes;i++){
432 for(
int j=0;j<numnodes;j++){
433 Ke->
values[i*numnodes+j] += D_scalar*h/(2.*vel)*dbasis[0*numnodes+i] *D[0]*dbasis[0*numnodes+j];
438 if(stabilization==2){
441 for(
int i=0;i<numnodes;i++){
442 for(
int j=0;j<numnodes;j++){
443 Ke->
values[i*numnodes+j]+=dt*gauss->
weight*Jdet*tau*(vx*dbasis[0*numnodes+i]+vy*dbasis[1*numnodes+i])*(vx*dbasis[0*numnodes+j]+vy*dbasis[1*numnodes+j]);
447 if(stabilization==5){
449 for(
int i=0;i<numnodes;i++){
450 for(
int j=0;j<numnodes;j++){
451 Ke->
values[i*numnodes+j]+=gauss->
weight*Jdet*tau*basis[j]*(vx*dbasis[0*numnodes+i]+vy*dbasis[1*numnodes+i]);
455 for(
int i=0;i<numnodes;i++){
456 for(
int j=0;j<numnodes;j++){
457 Ke->
values[i*numnodes+j]+=gauss->
weight*Jdet*tau*basis[j]*(basis[i]*dvxdx+basis[i]*dvydy);
462 for(
int i=0;i<numnodes;i++){
463 for(
int j=0;j<numnodes;j++){
464 Ke->
values[i*numnodes+j]+=dt*gauss->
weight*Jdet*tau*(vx*dbasis[0*numnodes+j]+vy*dbasis[1*numnodes+j])*(vx*dbasis[0*numnodes+i]+vy*dbasis[1*numnodes+i]);
468 for(
int i=0;i<numnodes;i++){
469 for(
int j=0;j<numnodes;j++){
470 Ke->
values[i*numnodes+j]+=dt*gauss->
weight*Jdet*tau*(vx*dbasis[0*numnodes+j]+vy*dbasis[1*numnodes+j])*(basis[i]*dvxdx+basis[i]*dvydy);
475 for(
int i=0;i<numnodes;i++){
476 for(
int j=0;j<numnodes;j++){
477 Ke->
values[i*numnodes+j]+=dt*gauss->
weight*Jdet*tau*(basis[j]*dvxdx+basis[j]*dvydy)*(vx*dbasis[0*numnodes+i]+vy*dbasis[1*numnodes+i]);
481 for(
int i=0;i<numnodes;i++){
482 for(
int j=0;j<numnodes;j++){
483 Ke->
values[i*numnodes+j]+=dt*gauss->
weight*Jdet*tau*(basis[j]*dvxdx+basis[j]*dvydy)*(basis[i]*dvxdx+basis[i]*dvydy);
490 xDelete<IssmDouble>(xyz_list);
491 xDelete<IssmDouble>(basis);
492 xDelete<IssmDouble>(dbasis);
◆ CreateKMatrixDG()
Definition at line 496 of file MasstransportAnalysis.cpp.
511 IssmDouble* basis = xNew<IssmDouble>(numnodes);
512 IssmDouble* dbasis = xNew<IssmDouble>(3*numnodes);
523 for(
int ig=gauss->
begin();ig<gauss->end();ig++){
533 D_scalar=gauss->
weight*Jdet;
534 for(
int i=0;i<numnodes;i++){
535 for(
int j=0;j<numnodes;j++){
536 Ke->
values[i*numnodes+j] += D_scalar*basis[i]*basis[j];
541 D_scalar = - dt*gauss->
weight*Jdet;
542 for(
int i=0;i<numnodes;i++){
543 for(
int j=0;j<numnodes;j++){
544 Ke->
values[i*numnodes+j] += D_scalar*(vx*dbasis[0*numnodes+i]*basis[j] + vy*dbasis[1*numnodes+i]*basis[j]);
550 xDelete<IssmDouble>(xyz_list);
551 xDelete<IssmDouble>(basis);
552 xDelete<IssmDouble>(dbasis);
◆ CreatePVector()
◆ CreatePVectorCG()
Definition at line 582 of file MasstransportAnalysis.cpp.
588 int stabilization,dim,domaintype;
589 int melt_style,point1;
614 IssmDouble* basis = xNew<IssmDouble>(numnodes);
615 IssmDouble* dbasis= xNew<IssmDouble>(dim*numnodes);
640 element->
GetGroundedPart(&point1,&fraction1,&fraction2,&mainlyfloating);
641 gauss = element->
NewGauss(point1,fraction1,fraction2,3);
648 for(
int ig=gauss->
begin();ig<gauss->end();ig++){
661 if (phi>0.999999999) mb=gmb;
662 else mb=(1-phi)*fmb+phi*gmb;
665 if(gllevelset>0.) mb=gmb;
669 if (phi<0.00000001) mb=fmb;
673 if (phi<0.99999999) mb=fmb;
678 for(
int i=0;i<numnodes;i++) pe->
values[i]+=Jdet*gauss->
weight*(thickness+dt*(ms-mb))*basis[i];
680 if(stabilization==5){
686 vel=sqrt(vx*vx+vy*vy)+1.e-8;
694 for(
int i=0;i<numnodes;i++){
695 pe->
values[i]+=Jdet*gauss->
weight*(thickness+dt*(ms-mb))*(tau*vx*dbasis[0*numnodes+i]+tau*vy*dbasis[1*numnodes+i]);
698 for(
int i=0;i<numnodes;i++){
699 pe->
values[i]+=Jdet*gauss->
weight*(thickness+dt*(ms-mb))*(tau*basis[i]*dvxdx+tau*basis[i]*dvydy);
706 xDelete<IssmDouble>(xyz_list);
707 xDelete<IssmDouble>(basis);
708 xDelete<IssmDouble>(dbasis);
◆ CreatePVectorDG()
Definition at line 712 of file MasstransportAnalysis.cpp.
718 int melt_style, point1;
730 IssmDouble* basis = xNew<IssmDouble>(numnodes);
746 element->
GetGroundedPart(&point1,&fraction1,&fraction2,&mainlyfloating);
747 gauss = element->
NewGauss(point1,fraction1,fraction2,3);
754 for(
int ig=gauss->
begin();ig<gauss->end();ig++){
767 if (phi>0.999999999) mb=gmb;
768 else mb=(1-phi)*fmb+phi*gmb;
771 if(gllevelset>0.) mb=gmb;
775 if (phi<0.00000001) mb=fmb;
779 if (phi<0.99999999) mb=fmb;
784 for(
int i=0;i<numnodes;i++) pe->
values[i]+=Jdet*gauss->
weight*(thickness+dt*(ms-mb))*basis[i];
788 xDelete<IssmDouble>(xyz_list);
789 xDelete<IssmDouble>(basis);
◆ GetSolutionFromInputs()
◆ GradientJ()
void MasstransportAnalysis::GradientJ |
( |
Vector< IssmDouble > * |
gradient, |
|
|
Element * |
element, |
|
|
int |
control_type, |
|
|
int |
control_index |
|
) |
| |
|
virtual |
◆ InputUpdateFromSolution()
void MasstransportAnalysis::InputUpdateFromSolution |
( |
IssmDouble * |
solution, |
|
|
Element * |
element |
|
) |
| |
|
virtual |
Implements Analysis.
Definition at line 799 of file MasstransportAnalysis.cpp.
809 IssmDouble* newthickness = xNew<IssmDouble>(numnodes);
810 IssmDouble* thicknessresidual = xNew<IssmDouble>(numnodes);
814 for(
int i=0;i<numnodes;i++){
815 newthickness[i]=solution[doflist[i]];
816 thicknessresidual[i]=0.;
818 if(xIsNan<IssmDouble>(newthickness[i]))
_error_(
"NaN found in solution vector");
819 if(xIsInf<IssmDouble>(newthickness[i]))
_error_(
"Inf found in solution vector");
820 if(newthickness[i]<minthickness){
821 thicknessresidual[i]=minthickness-newthickness[i];
822 newthickness[i]=minthickness;
828 xDelete<int>(doflist);
829 xDelete<IssmDouble>(newthickness);
830 xDelete<IssmDouble>(thicknessresidual);
843 newthickness = xNew<IssmDouble>(numvertices);
844 IssmDouble* oldthickness = xNew<IssmDouble>(numvertices);
845 IssmDouble* cumdeltathickness = xNew<IssmDouble>(numvertices);
846 IssmDouble* deltathickness = xNew<IssmDouble>(numvertices);
847 IssmDouble* newbase = xNew<IssmDouble>(numvertices);
848 IssmDouble* bed = xNew<IssmDouble>(numvertices);
849 IssmDouble* newsurface = xNew<IssmDouble>(numvertices);
850 IssmDouble* oldbase = xNew<IssmDouble>(numvertices);
851 IssmDouble* oldsurface = xNew<IssmDouble>(numvertices);
852 IssmDouble* phi = xNew<IssmDouble>(numvertices);
853 IssmDouble* sealevel = xNew<IssmDouble>(numvertices);
865 bool isgroundingline;
870 for(
int i=0;i<numvertices;i++){
871 cumdeltathickness[i] += newthickness[i]-oldthickness[i];
872 deltathickness[i] = newthickness[i]-oldthickness[i];
881 for(
int i=0;i<numvertices;i++) {
884 newsurface[i] = bed[i]+newthickness[i];
888 newsurface[i] = oldbase[i]+newthickness[i];
889 newbase[i] = oldbase[i];
894 newsurface[i] = newthickness[i]*(1.-rho_ice/rho_water)+sealevel[i];
895 newbase[i] = newthickness[i]*(-rho_ice/rho_water)+sealevel[i];
898 newsurface[i] = oldsurface[i]+(1.0-rho_ice/rho_water)*(newthickness[i]-oldthickness[i])+sealevel[i];
899 newbase[i] = oldbase[i]-rho_ice/rho_water*(newthickness[i]-oldthickness[i])+sealevel[i];
901 else _error_(
"Hydrostatic adjustment " << hydroadjustment <<
" (" <<
EnumToStringx(hydroadjustment) <<
") not supported yet");
912 xDelete<IssmDouble>(newthickness);
913 xDelete<IssmDouble>(cumdeltathickness);
914 xDelete<IssmDouble>(newbase);
915 xDelete<IssmDouble>(newsurface);
916 xDelete<IssmDouble>(oldthickness);
917 xDelete<IssmDouble>(deltathickness);
918 xDelete<IssmDouble>(oldbase);
919 xDelete<IssmDouble>(oldsurface);
920 xDelete<IssmDouble>(phi);
921 xDelete<IssmDouble>(sealevel);
922 xDelete<IssmDouble>(bed);
923 xDelete<int>(doflist);
◆ UpdateConstraints()
void MasstransportAnalysis::UpdateConstraints |
( |
FemModel * |
femmodel | ) |
|
|
virtual |
◆ CreateFctKMatrix()
Definition at line 931 of file MasstransportAnalysis.cpp.
948 IssmDouble* basis = xNew<IssmDouble>(numnodes);
949 IssmDouble* dbasis = xNew<IssmDouble>(dim*numnodes);
950 IssmDouble* D = xNewZeroInit<IssmDouble>(dim*dim);
959 for(
int ig=gauss->
begin();ig<gauss->end();ig++){
972 D_scalar = gauss->
weight*Jdet;
973 for(
int i=0;i<numnodes;i++){
974 for(
int j=0;j<numnodes;j++){
976 Ke->
values[i*numnodes+j] += -D_scalar*basis[i]*basis[j]*(dvxdx+dvydy);
978 Ke->
values[i*numnodes+j] += -D_scalar*(vx*dbasis[0*numnodes+j]*basis[i] + vy*dbasis[1*numnodes+j]*basis[i]);
984 xDelete<IssmDouble>(xyz_list);
985 xDelete<IssmDouble>(D);
986 xDelete<IssmDouble>(basis);
987 xDelete<IssmDouble>(dbasis);
◆ CreateMassMatrix()
Definition at line 991 of file MasstransportAnalysis.cpp.
1005 IssmDouble* basis = xNew<IssmDouble>(numnodes);
1012 for(
int ig=gauss->
begin();ig<gauss->end();ig++){
1026 xDelete<IssmDouble>(xyz_list);
1027 xDelete<IssmDouble>(basis);
◆ CreateFctPVector()
Definition at line 1031 of file MasstransportAnalysis.cpp.
1037 int stabilization,dim,domaintype;
1038 int melt_style,point1;
1039 bool mainlyfloating;
1047 Gauss* gauss = NULL;
1063 IssmDouble* basis = xNew<IssmDouble>(numnodes);
1064 IssmDouble* dbasis= xNew<IssmDouble>(dim*numnodes);
1080 element->
GetGroundedPart(&point1,&fraction1,&fraction2,&mainlyfloating);
1081 gauss = element->
NewGauss(point1,fraction1,fraction2,3);
1088 for(
int ig=gauss->
begin();ig<gauss->end();ig++){
1100 if (phi>0.999999999) mb=gmb;
1101 else mb=(1-phi)*fmb+phi*gmb;
1104 if(gllevelset>0.) mb=gmb;
1108 if (phi<0.00000001) mb=fmb;
1112 if (phi<0.99999999) mb=fmb;
1117 for(
int i=0;i<numnodes;i++) pe->
values[i]+=Jdet*gauss->
weight*(ms-mb)*basis[i];
1122 xDelete<IssmDouble>(xyz_list);
1123 xDelete<IssmDouble>(basis);
1124 xDelete<IssmDouble>(dbasis);
◆ FctKMatrix()
◆ FctPVector()
◆ MassMatrix()
The documentation for this class was generated from the following files:
@ BasalforcingsPicoBasinIdEnum
@ BasalforcingsOceanTempEnum
virtual int GetElementType(void)=0
ElementMatrix * CreateFctKMatrix(Element *element)
virtual int GetNumberOfNodes(void)=0
@ MasstransportAnalysisEnum
void FindParam(bool *pvalue, int paramenum)
int AddObject(Object *object)
@ TimesteppingTimeStepEnum
virtual Element * SpawnBasalElement(void)=0
int TripleMultiply(IssmDouble *a, int nrowa, int ncola, int itrna, IssmDouble *b, int nrowb, int ncolb, int itrnb, IssmDouble *c, int nrowc, int ncolc, int itrnc, IssmDouble *d, int iaddd)
@ MismipFloatingMeltRateEnum
void AllocateSystemMatricesx(Matrix< IssmDouble > **pKff, Matrix< IssmDouble > **pKfs, Vector< IssmDouble > **pdf, Vector< IssmDouble > **ppf, FemModel *femmodel)
@ MasstransportStabilizationEnum
@ MasstransportRequestedOutputsEnum
@ BasalforcingsUpperwaterElevationEnum
@ TransientIsgroundinglineEnum
@ BasalforcingsDeepwaterMeltingRateEnum
void AddObject(Param *newparam)
ElementVector * CreateFctPVector(Element *element)
virtual Input2 * GetInput2(int inputenum)=0
void DeleteMaterials(void)
virtual void NodalFunctions(IssmDouble *basis, Gauss *gauss)=0
@ SealevelriseCumDeltathicknessEnum
ElementVector * NewElementVector(int approximation_enum=NoneApproximationEnum)
void DeleteData(int num,...)
void ActivatePenaltyMethod(int in_analysis)
@ BasalforcingsPerturbationMeltingRateEnum
Param * CopyConstantObject(const char *constant_name, int param_enum)
@ MasstransportSpcthicknessEnum
@ BasalforcingsGroundediceMeltingRateEnum
virtual void GetGroundedPart(int *point1, IssmDouble *fraction1, IssmDouble *fraction2, bool *mainlyfloating)=0
@ NoMeltOnPartiallyFloatingEnum
@ FullMeltOnPartiallyFloatingEnum
virtual Gauss * NewGauss(void)=0
@ BasalforcingsOceanSalinityEnum
const char * EnumToStringx(int enum_in)
@ LinearFloatingMeltRateEnum
@ BeckmannGoosseFloatingMeltRateEnum
@ MantlePlumeGeothermalFluxEnum
@ BasalforcingsFloatingiceMeltingRateEnum
void FindConstant(bool *pvalue, const char *constant_name)
@ GroundinglineMeltInterpolationEnum
void GetVerticesCoordinates(IssmDouble **xyz_list)
@ BasalforcingsPicoOverturningCoeffEnum
@ SpatialLinearFloatingMeltRateEnum
void FetchData(bool *pboolean, const char *data_name)
@ MaterialsRhoSeawaterEnum
@ BasalforcingsIsmip6Enum
ElementMatrix * CreateKMatrixCG(Element *element)
@ MasstransportMinThicknessEnum
void GetDofListLocal(int **pdoflist, int approximation_enum, int setenum)
@ BasalforcingsIsmip6TfEnum
virtual void AddBasalInput2(int input_enum, IssmDouble *values, int interpolation_enum)
void CreateNodes(Nodes *nodes, IoModel *iomodel, bool isamr=false)
@ BasalforcingsIsmip6BasinIdEnum
ElementVector * CreatePVectorCG(Element *element)
@ MasstransportIsfreesurfaceEnum
void GetSolutionFromInputsOneDof(Vector< IssmDouble > *solution, int solutionenum)
ElementVector * CreatePVectorDG(Element *element)
@ BasalforcingsDeepwaterElevationEnum
#define _error_(StreamArgs)
void SetActiveNodesLSMx(FemModel *femmodel)
virtual int begin(void)=0
Object * GetObjectByOffset(int offset)
virtual IssmDouble CharacteristicLength(void)=0
@ MeshVertexonsurfaceEnum
@ SurfaceloadIceThicknessChangeEnum
virtual IssmDouble GetGroundedPortion(IssmDouble *xyz_list)=0
virtual void JacobianDeterminant(IssmDouble *Jdet, IssmDouble *xyz_list, Gauss *gauss)=0
virtual void GaussPoint(int ig)=0
virtual void Update(Inputs2 *inputs2, int index, IoModel *iomodel, int analysis_counter, int analysis_type, int finite_element)=0
ElementMatrix * CreateKMatrixDG(Element *element)
@ BasalforcingsIsmip6MeltAnomalyEnum
virtual int FiniteElement(void)=0
void AddToGlobal(Vector< IssmDouble > *pf)
void FetchDataToInput(Inputs2 *inputs2, Elements *elements, const char *vector_name, int input_enum)
void IoModelToConstraintsx(Constraints *constraints, IoModel *iomodel, const char *spc_name, int analysis_type, int finite_element, int dof)
void AddToGlobal(Matrix< IssmDouble > *Kff, Matrix< IssmDouble > *Kfs)
@ SealevelriseCumDeltathicknessOldEnum
void CreateFaces(IoModel *iomodel)
void DatasetInputAdd(int enum_type, IssmDouble *vector, Inputs2 *inputs2, IoModel *iomodel, int M, int N, int vector_type, int vector_enum, int code, int input_enum)
@ MasstransportNumRequestedOutputsEnum
virtual int GetNumberOfVertices(void)=0
void GetInputListOnVertices(IssmDouble *pvalue, int enumtype)
@ MasstransportPenaltyFactorEnum
virtual void NodalFunctionsDerivatives(IssmDouble *dbasis, IssmDouble *xyz_list, Gauss *gauss)=0
ElementMatrix * CreateMassMatrix(Element *element)
@ MasstransportHydrostaticAdjustmentEnum
ElementMatrix * NewElementMatrix(int approximation_enum=NoneApproximationEnum)