source: issm/oecreview/Archive/16554-17801/ISSM-17393-17394.diff

Last change on this file was 17802, checked in by Mathieu Morlighem, 11 years ago

Added archives
File size: 7.5 KB
RevLine 
[17802]1Index: ../trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp
2===================================================================
3--- ../trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp (revision 17393)
4+++ ../trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp (revision 17394)
5@@ -6,7 +6,7 @@
6 #include "../solutionsequences/solutionsequences.h"
7 #include "../cores/cores.h"
8
9-//#define FSANALYTICAL 101
10+//#define FSANALYTICAL 4
11
12 /*Model processing*/
13 int StressbalanceAnalysis::DofsPerNode(int** pdoftype,int meshtype,int approximation){/*{{{*/
14Index: ../trunk-jpl/src/c/shared/Enum/EnumDefinitions.h
15===================================================================
16--- ../trunk-jpl/src/c/shared/Enum/EnumDefinitions.h (revision 17393)
17+++ ../trunk-jpl/src/c/shared/Enum/EnumDefinitions.h (revision 17394)
18@@ -532,6 +532,7 @@
19 StressTensoryyEnum,
20 StressTensoryzEnum,
21 StressTensorzzEnum,
22+ DivergenceEnum,
23 GiaCrossSectionShapeEnum,
24 GiadWdtEnum,
25 GiaWEnum,
26Index: ../trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp
27===================================================================
28--- ../trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp (revision 17393)
29+++ ../trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp (revision 17394)
30@@ -521,6 +521,7 @@
31 case StressTensoryyEnum : return "StressTensoryy";
32 case StressTensoryzEnum : return "StressTensoryz";
33 case StressTensorzzEnum : return "StressTensorzz";
34+ case DivergenceEnum : return "Divergence";
35 case GiaCrossSectionShapeEnum : return "GiaCrossSectionShape";
36 case GiadWdtEnum : return "GiadWdt";
37 case GiaWEnum : return "GiaW";
38Index: ../trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp
39===================================================================
40--- ../trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp (revision 17393)
41+++ ../trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp (revision 17394)
42@@ -533,6 +533,7 @@
43 else if (strcmp(name,"StressTensoryy")==0) return StressTensoryyEnum;
44 else if (strcmp(name,"StressTensoryz")==0) return StressTensoryzEnum;
45 else if (strcmp(name,"StressTensorzz")==0) return StressTensorzzEnum;
46+ else if (strcmp(name,"Divergence")==0) return DivergenceEnum;
47 else if (strcmp(name,"GiaCrossSectionShape")==0) return GiaCrossSectionShapeEnum;
48 else if (strcmp(name,"GiadWdt")==0) return GiadWdtEnum;
49 else if (strcmp(name,"GiaW")==0) return GiaWEnum;
50@@ -627,11 +628,11 @@
51 else if (strcmp(name,"QmuInName")==0) return QmuInNameEnum;
52 else if (strcmp(name,"QmuOutName")==0) return QmuOutNameEnum;
53 else if (strcmp(name,"Regular")==0) return RegularEnum;
54- else if (strcmp(name,"Scaled")==0) return ScaledEnum;
55 else stage=6;
56 }
57 if(stage==6){
58- if (strcmp(name,"Separate")==0) return SeparateEnum;
59+ if (strcmp(name,"Scaled")==0) return ScaledEnum;
60+ else if (strcmp(name,"Separate")==0) return SeparateEnum;
61 else if (strcmp(name,"Sset")==0) return SsetEnum;
62 else if (strcmp(name,"Verbose")==0) return VerboseEnum;
63 else if (strcmp(name,"TriangleInterp")==0) return TriangleInterpEnum;
64Index: ../trunk-jpl/src/c/classes/FemModel.cpp
65===================================================================
66--- ../trunk-jpl/src/c/classes/FemModel.cpp (revision 17393)
67+++ ../trunk-jpl/src/c/classes/FemModel.cpp (revision 17394)
68@@ -431,6 +431,7 @@
69
70 switch (response_descriptor_enum){
71
72+ case DivergenceEnum: this->Divergencex(responses); break;
73 case IceVolumeEnum: this->IceVolumex(responses); break;
74 case IceVolumeAboveFloatationEnum: this->IceVolumeAboveFloatationx(responses); break;
75 case MinVelEnum: this->MinVelx(responses); break;
76@@ -503,6 +504,7 @@
77 switch(output_enum){
78
79 /*Scalar output*/
80+ case DivergenceEnum: this->Divergencex(&double_result); break;
81 case IceVolumeEnum: this->IceVolumex(&double_result); break;
82 case IceVolumeAboveFloatationEnum: this->IceVolumeAboveFloatationx(&double_result); break;
83 case MinVelEnum: this->MinVelx(&double_result); break;
84@@ -1041,6 +1043,22 @@
85 *pSmb=total_smb;
86
87 }/*}}}*/
88+void FemModel::Divergencex(IssmDouble* pdiv){/*{{{*/
89+
90+ IssmDouble local_divergence=0;
91+ IssmDouble total_divergence;
92+
93+ for(int i=0;i<this->elements->Size();i++){
94+ Element* element=dynamic_cast<Element*>(this->elements->GetObjectByOffset(i));
95+ local_divergence+=element->Divergence();
96+ }
97+ ISSM_MPI_Reduce(&local_divergence,&total_divergence,1,ISSM_MPI_DOUBLE,ISSM_MPI_SUM,0,IssmComm::GetComm() );
98+ ISSM_MPI_Bcast(&total_divergence,1,ISSM_MPI_DOUBLE,0,IssmComm::GetComm());
99+
100+ /*Assign output pointers: */
101+ *pdiv=total_divergence;
102+
103+}/*}}}*/
104 void FemModel::IceVolumex(IssmDouble* pV){/*{{{*/
105
106 IssmDouble local_ice_volume = 0;
107Index: ../trunk-jpl/src/c/classes/Elements/Element.h
108===================================================================
109--- ../trunk-jpl/src/c/classes/Elements/Element.h (revision 17393)
110+++ ../trunk-jpl/src/c/classes/Elements/Element.h (revision 17394)
111@@ -55,6 +55,7 @@
112 /* bool AnyActive(void); */
113 void CoordinateSystemTransform(IssmDouble** ptransform,Node** nodes,int numnodes,int* cs_array);
114 void DeleteMaterials(void);
115+ IssmDouble Divergence(void);
116 void FindParam(bool* pvalue,int paramenum);
117 void FindParam(int* pvalue,int paramenum);
118 void FindParam(IssmDouble* pvalue,int paramenum);
119Index: ../trunk-jpl/src/c/classes/Elements/Element.cpp
120===================================================================
121--- ../trunk-jpl/src/c/classes/Elements/Element.cpp (revision 17393)
122+++ ../trunk-jpl/src/c/classes/Elements/Element.cpp (revision 17394)
123@@ -129,6 +129,38 @@
124 void Element::DeleteMaterials(void){/*{{{*/
125 delete this->material;
126 }/*}}}*/
127+IssmDouble Element::Divergence(void){/*{{{*/
128+ /*Compute element divergence*/
129+
130+ /*Intermediaries*/
131+ IssmDouble Jdet;
132+ IssmDouble divergence=0.;
133+ IssmDouble dvx[3],dvy[3],dvz[3];
134+ IssmDouble *xyz_list = NULL;
135+
136+ /*Get inputs and parameters*/
137+ Input* vx_input = this->GetInput(VxEnum); _assert_(vx_input);
138+ Input* vy_input = this->GetInput(VyEnum); _assert_(vy_input);
139+ Input* vz_input = this->GetInput(VzEnum); _assert_(vz_input);
140+ this->GetVerticesCoordinates(&xyz_list);
141+
142+ Gauss* gauss=this->NewGaussBase(5);
143+ for(int ig=gauss->begin();ig<gauss->end();ig++){
144+ gauss->GaussPoint(ig);
145+ this->JacobianDeterminant(&Jdet,xyz_list,gauss);
146+
147+ /*Get strain rate assuming that epsilon has been allocated*/
148+ vx_input->GetInputDerivativeValue(&dvx[0],xyz_list,gauss);
149+ vy_input->GetInputDerivativeValue(&dvy[0],xyz_list,gauss);
150+ vz_input->GetInputDerivativeValue(&dvz[0],xyz_list,gauss);
151+
152+ divergence += (dvx[0]+dvy[1]+dvz[2])*gauss->weight*Jdet;
153+ }
154+
155+ /*Clean up and return*/
156+ delete gauss;
157+ return divergence;
158+}/*}}}*/
159 void Element::FindParam(bool* pvalue,int paramenum){/*{{{*/
160 this->parameters->FindParam(pvalue,paramenum);
161 }/*}}}*/
162Index: ../trunk-jpl/src/c/classes/FemModel.h
163===================================================================
164--- ../trunk-jpl/src/c/classes/FemModel.h (revision 17393)
165+++ ../trunk-jpl/src/c/classes/FemModel.h (revision 17394)
166@@ -71,6 +71,7 @@
167 void MinVyx(IssmDouble* presponse);
168 void MinVzx(IssmDouble* presponse);
169 void TotalSmbx(IssmDouble* pSmb);
170+ void Divergencex(IssmDouble* pdiv);
171 void IceVolumex(IssmDouble* pV);
172 void IceVolumeAboveFloatationx(IssmDouble* pV);
173 void ElementResponsex(IssmDouble* presponse,int response_enum);
Note: See TracBrowser for help on using the repository browser.