Changeset 26223

Timestamp:

04/29/21 13:31:33 (4 years ago)

Author:

Mathieu Morlighem

Message:

CHG: trying to fix Eric's commit

Location:

issm/trunk-jpl/src/c

Files:

• analyses/DamageEvolutionAnalysis.cpp (modified) (2 diffs)
• analyses/DamageEvolutionAnalysis.h (modified) (1 diff)
• classes/Elements/Tria.cpp (modified) (1 diff)
• classes/SealevelGeometry.cpp (modified) (2 diffs)
• classes/SealevelGeometry.h (modified) (2 diffs)
• solutionsequences/solutionsequence_fct.cpp (modified) (1 diff)

issm/trunk-jpl/src/c/analyses/DamageEvolutionAnalysis.cpp

@@ -809 +809 @@
         return Me;
 }/*}}}*/
+ElementVector* DamageEvolutionAnalysis::CreateFctPVector(Element* element){/*{{{*/
+
+        return this->CreatePVector(element);
+
+}/*}}}*/
 void           DamageEvolutionAnalysis::FctKMatrix(Matrix<IssmDouble>** pKff,Matrix<IssmDouble>** pKfs,FemModel* femmodel){/*{{{*/
 
@@ -837 +842 @@
         }
 }/*}}}*/
+void           DamageEvolutionAnalysis::FctPVector(Vector<IssmDouble>** ppf,FemModel* femmodel){/*{{{*/
+
+        /*Output*/
+        Vector<IssmDouble>* pf = NULL;
+
+        /*Initialize P vector*/
+        AllocateSystemMatricesx(NULL,NULL,NULL,&pf,femmodel);
+
+        /*Create and assemble matrix*/
+        for(Object* & object : femmodel->elements->objects){
+                Element*       element = xDynamicCast<Element*>(object);
+                ElementVector* pe      = this->CreateFctPVector(element);
+                if(pe) pe->AddToGlobal(pf);
+                delete pe;
+        }
+        pf->Assemble();
+
+        /*Assign output pointer*/
+        *ppf=pf;
+}/*}}}*/
 void           DamageEvolutionAnalysis::MassMatrix(Matrix<IssmDouble>** pMff,FemModel* femmodel){/*{{{*/
 

issm/trunk-jpl/src/c/analyses/DamageEvolutionAnalysis.h

@@ -39 +39 @@
                 ElementMatrix* CreateFctKMatrix(Element* element);
                 ElementMatrix* CreateMassMatrix(Element* element);
+                ElementVector* CreateFctPVector(Element* element);
                 void           FctKMatrix(Matrix<IssmDouble>** pKff,Matrix<IssmDouble>** pKfs,FemModel* femmodel);
                 void           MassMatrix(Matrix<IssmDouble>** pMff,FemModel* femmodel);
+                void                            FctPVector(Vector<IssmDouble>** ppf,FemModel* femmodel);
 };
 #endif

issm/trunk-jpl/src/c/classes/Elements/Tria.cpp

@@ -1816 +1816 @@
         }
         //normalize to phi. 
-        if(total_weight)for(int i=0;i<NUMVERTICES;i++)weights[i]/=total_weight/phi; 
+        if(total_weight>0.) for(int i=0;i<NUMVERTICES;i++)weights[i]/=total_weight/phi; 
         else for(int i=0;i<NUMVERTICES;i++)weights[i]=0;
 

issm/trunk-jpl/src/c/classes/SealevelGeometry.cpp

@@ -20 +20 @@
         localnel=localnelin;
         for(int i=0;i<SLGEOM_NUMLOADS;i++){
-                for (int j=0;j<MAXVERTICES;j++){
+                for (int j=0;j<SLMAXVERTICES;j++){
                         LoadWeigths[i][j]=xNewZeroInit<IssmDouble>(localnel);
                 }
@@ -43 +43 @@
 SealevelGeometry::~SealevelGeometry(){ /*{{{*/
         for(int i=0;i<SLGEOM_NUMLOADS;i++){
-                for (int j=0;j<MAXVERTICES;j++){
+                for (int j=0;j<SLMAXVERTICES;j++){
                         xDelete<IssmDouble>(LoadWeigths[i][j]);
                 }

issm/trunk-jpl/src/c/classes/SealevelGeometry.h

@@ -11 +11 @@
 #define SLGEOM_ICE 1 
 #define SLGEOM_WATER 2
-#define MAXVERTICES 3
+#define SLMAXVERTICES 3
 
 #include "../toolkits/toolkits.h"
@@ -20 +20 @@
 
                 int         localnel;
-                IssmDouble* LoadWeigths[SLGEOM_NUMLOADS][MAXVERTICES];
+                IssmDouble* LoadWeigths[SLGEOM_NUMLOADS][SLMAXVERTICES];
                 IssmDouble* LoadArea[SLGEOM_NUMLOADS];
                 Vector<IssmDouble>* vlatbarycentre[SLGEOM_NUMLOADS];

issm/trunk-jpl/src/c/solutionsequences/solutionsequence_fct.cpp

@@ -418 +418 @@
                         danalysis->MassMatrix(&Mc,femmodel);
                         danalysis->FctKMatrix(&K,NULL,femmodel);
+                        danalysis->FctPVector(&p,femmodel);
                         break;
                 default: _error_("analysis type " << EnumToStringx(analysis_type) << " not supported for FCT\n");