source: issm/oecreview/Archive/15392-16133/ISSM-15813-15814.diff@ 16134

Last change on this file since 16134 was 16134, checked in by Mathieu Morlighem, 12 years ago

Added Archive/15392-16133
File size: 5.7 KB
RevLine 
[16134]1Index: ../trunk-jpl/src/c/classes/Elements/Tria.cpp
2===================================================================
3--- ../trunk-jpl/src/c/classes/Elements/Tria.cpp (revision 15813)
4+++ ../trunk-jpl/src/c/classes/Elements/Tria.cpp (revision 15814)
5@@ -817,16 +817,16 @@
6 bool floating=true;
7 int point;
8 const IssmPDouble epsilon= 1.e-15;
9- IssmDouble gl[3];
10+ IssmDouble gl[NUMVERTICES];
11 IssmDouble f1,f2;
12
13 /*Recover parameters and values*/
14 GetInputListOnVertices(&gl[0],GLlevelsetEnum);
15
16 /*Be sure that values are not zero*/
17- if(gl[0]==0) gl[0]=gl[0]+epsilon;
18- if(gl[1]==0) gl[1]=gl[1]+epsilon;
19- if(gl[2]==0) gl[2]=gl[2]+epsilon;
20+ if(gl[0]==0.) gl[0]=gl[0]+epsilon;
21+ if(gl[1]==0.) gl[1]=gl[1]+epsilon;
22+ if(gl[2]==0.) gl[2]=gl[2]+epsilon;
23
24 /*Check that not all nodes are grounded or floating*/
25 if(gl[0]>0 && gl[1]>0 && gl[2]>0){ // All grounded
26@@ -868,19 +868,19 @@
27 IssmDouble Tria::GetGroundedPortion(IssmDouble* xyz_list){
28 /*Computeportion of the element that is grounded*/
29
30- bool mainlyfloating = true;
31- const IssmPDouble epsilon= 1.e-15;
32+ bool mainlyfloating = true;
33+ const IssmPDouble epsilon = 1.e-15;
34 IssmDouble phi,s1,s2,area_init,area_grounded;
35- IssmDouble gl[3];
36+ IssmDouble gl[NUMVERTICES];
37 IssmDouble xyz_bis[3][3];
38
39 /*Recover parameters and values*/
40 GetInputListOnVertices(&gl[0],GLlevelsetEnum);
41
42 /*Be sure that values are not zero*/
43- if(gl[0]==0) gl[0]=gl[0]+epsilon;
44- if(gl[1]==0) gl[1]=gl[1]+epsilon;
45- if(gl[2]==0) gl[2]=gl[2]+epsilon;
46+ if(gl[0]==0.) gl[0]=gl[0]+epsilon;
47+ if(gl[1]==0.) gl[1]=gl[1]+epsilon;
48+ if(gl[2]==0.) gl[2]=gl[2]+epsilon;
49
50 /*Check that not all nodes are grounded or floating*/
51 if(gl[0]>0 && gl[1]>0 && gl[2]>0){ // All grounded
52Index: ../trunk-jpl/src/c/classes/Elements/Penta.cpp
53===================================================================
54--- ../trunk-jpl/src/c/classes/Elements/Penta.cpp (revision 15813)
55+++ ../trunk-jpl/src/c/classes/Elements/Penta.cpp (revision 15814)
56@@ -1023,16 +1023,16 @@
57 bool mainlyfloating = true;
58 const IssmPDouble epsilon= 1.e-15;
59 IssmDouble phi,s1,s2,area_init,area_grounded;
60- IssmDouble gl[3];
61- IssmDouble xyz_bis[3][3];
62+ IssmDouble gl[NUMVERTICES];
63+ IssmDouble xyz_bis[NUMVERTICES2D][3];
64
65 /*Recover parameters and values*/
66 GetInputListOnVertices(&gl[0],GLlevelsetEnum);
67
68 /*Be sure that values are not zero*/
69- if(gl[0]==0) gl[0]=gl[0]+epsilon;
70- if(gl[1]==0) gl[1]=gl[1]+epsilon;
71- if(gl[2]==0) gl[2]=gl[2]+epsilon;
72+ if(gl[0]==0.) gl[0]=gl[0]+epsilon;
73+ if(gl[1]==0.) gl[1]=gl[1]+epsilon;
74+ if(gl[2]==0.) gl[2]=gl[2]+epsilon;
75
76 /*Check that not all nodes are grounded or floating*/
77 if(gl[0]>0 && gl[1]>0 && gl[2]>0){ // All grounded
78@@ -1047,23 +1047,23 @@
79
80 if(gl[0]*gl[1]>0){ //Nodes 0 and 1 are similar, so points must be found on segment 0-2 and 1-2
81 /*Coordinates of point 2: same as initial point 2*/
82- xyz_bis[2][0]=*(xyz_list+3*2+0);
83- xyz_bis[2][1]=*(xyz_list+3*2+1);
84- xyz_bis[2][2]=*(xyz_list+3*2+2);
85+ xyz_bis[2][0]=xyz_list[3*2+0];
86+ xyz_bis[2][1]=xyz_list[3*2+1];
87+ xyz_bis[2][2]=xyz_list[3*2+2];
88
89 /*Portion of the segments*/
90 s1=gl[2]/(gl[2]-gl[1]);
91 s2=gl[2]/(gl[2]-gl[0]);
92
93 /*New point 1*/
94- xyz_bis[1][0]=*(xyz_list+3*2+0)+s1*(*(xyz_list+3*1+0)-*(xyz_list+3*2+0));
95- xyz_bis[1][1]=*(xyz_list+3*2+1)+s1*(*(xyz_list+3*1+1)-*(xyz_list+3*2+1));
96- xyz_bis[1][2]=*(xyz_list+3*2+2)+s1*(*(xyz_list+3*1+2)-*(xyz_list+3*2+2));
97+ xyz_bis[1][0]=xyz_list[3*2+0]+s1*(xyz_list[3*1+0]-xyz_list[3*2+0]);
98+ xyz_bis[1][1]=xyz_list[3*2+1]+s1*(xyz_list[3*1+1]-xyz_list[3*2+1]);
99+ xyz_bis[1][2]=xyz_list[3*2+2]+s1*(xyz_list[3*1+2]-xyz_list[3*2+2]);
100
101 /*New point 0*/
102- xyz_bis[0][0]=*(xyz_list+3*2+0)+s2*(*(xyz_list+3*0+0)-*(xyz_list+3*2+0));
103- xyz_bis[0][1]=*(xyz_list+3*2+1)+s2*(*(xyz_list+3*0+1)-*(xyz_list+3*2+1));
104- xyz_bis[0][2]=*(xyz_list+3*2+2)+s2*(*(xyz_list+3*0+2)-*(xyz_list+3*2+2));
105+ xyz_bis[0][0]=xyz_list[3*2+0]+s2*(xyz_list[3*0+0]-xyz_list[3*2+0]);
106+ xyz_bis[0][1]=xyz_list[3*2+1]+s2*(xyz_list[3*0+1]-xyz_list[3*2+1]);
107+ xyz_bis[0][2]=xyz_list[3*2+2]+s2*(xyz_list[3*0+2]-xyz_list[3*2+2]);
108 }
109 else if(gl[1]*gl[2]>0){ //Nodes 1 and 2 are similar, so points must be found on segment 0-1 and 0-2
110 /*Coordinates of point 0: same as initial point 2*/
111@@ -1107,13 +1107,13 @@
112 }
113
114 /*Compute fraction of grounded element*/
115- GetJacobianDeterminant(&area_init, xyz_list,NULL);
116- GetJacobianDeterminant(&area_grounded, &xyz_bis[0][0],NULL);
117+ GetTriaJacobianDeterminant(&area_init, xyz_list,NULL);
118+ GetTriaJacobianDeterminant(&area_grounded, &xyz_bis[0][0],NULL);
119 if(mainlyfloating==true) area_grounded=area_init-area_grounded;
120 phi=area_grounded/area_init;
121 }
122
123- if(phi>1 || phi<0) _error_("Error. Problem with portion of grounded element: value should be between 0 and 1");
124+ if(phi>1. || phi<0.) _error_("Error. Problem with portion of grounded element: value should be between 0 and 1");
125
126 return phi;
127 }
128@@ -7659,7 +7659,7 @@
129 int i,j;
130 int analysis_type,migration_style;
131 IssmDouble xyz_list[NUMVERTICES][3];
132- IssmDouble xyz_list_tria[NUMVERTICES2D][3]={0.};
133+ IssmDouble xyz_list_tria[NUMVERTICES2D][3];
134 IssmDouble alpha2,Jdet;
135 IssmDouble phi=1.0;
136 IssmDouble DL_scalar;
137@@ -7692,7 +7692,7 @@
138 friction=new Friction("2d",inputs,matpar,analysis_type);
139
140 /*Recover portion of element that is grounded*/
141- if(migration_style==SubelementMigrationEnum) phi=this->GetGroundedPortion(&xyz_list[0][0]);
142+ if(migration_style==SubelementMigrationEnum) phi=this->GetGroundedPortion(&xyz_list_tria[0][0]);
143
144 /* Start looping on the number of gaussian points: */
145 gauss=new GaussPenta(0,1,2,2);
Note: See TracBrowser for help on using the repository browser.