Changeset 12976

Timestamp:

08/10/12 10:41:52 (13 years ago)

Author:

Mathieu Morlighem

Message:

normalize velocity

File:

• issm/trunk-jpl/src/c/classes/objects/Elements/Tria.cpp (modified) (6 diffs)

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

@@ -4305 +4305 @@
         IssmDouble  xyz_list[NUMVERTICES][3];
         IssmDouble  dp[NDOF2];
-        IssmDouble  vx,vy;
+        IssmDouble  vx,vy,vel;
         GaussTria  *gauss                    = NULL;
 
@@ -4334 +4334 @@
                 vx_input->GetInputValue(&vx,gauss);
                 vy_input->GetInputValue(&vy,gauss);
+                vel = sqrt(vx*vx+vy*vy);
+                vx  = vx/(vel+1.e-9);
+                vy  = vy/(vel+1.e-9);
 
                 /*J = 1/2 ( vx*dH/dx + vy*dH/dy )^2 */
@@ -4354 +4357 @@
         IssmDouble  xyz_list[NUMVERTICES][3];
         IssmDouble  dp[NDOF2];
-        IssmDouble  vx,vy;
+        IssmDouble  vx,vy,vel;
         GaussTria  *gauss                    = NULL;
 
@@ -4383 +4386 @@
                 vx_input->GetInputValue(&vx,gauss);
                 vy_input->GetInputValue(&vy,gauss);
+                vel = sqrt(vx*vx+vy*vy);
+                vx  = vx/(vel+1.e-9);
+                vy  = vy/(vel+1.e-9);
 
                 /*J = 1/2 ( -vy*dH/dx + vx*dH/dy )^2 */
@@ -4453 +4459 @@
         IssmDouble  dbasis[NDOF2][NUMVERTICES];
         IssmDouble  dH[2];
-        IssmDouble  v[2];
+        IssmDouble  vx,vy,vel;
         GaussTria *gauss     = NULL;
         int       *responses = NULL;
@@ -4498 +4504 @@
                         case ThicknessAlongGradientEnum:
                                 weights_input->GetInputValue(&weight, gauss,resp);
-                                vx_input->GetInputValue(&v[0],gauss);
-                                vy_input->GetInputValue(&v[1],gauss);
-                                for(i=0;i<numdof;i++) pe->values[i]+= - weight*(dH[0]*v[0]+dH[1]*v[1])*(dbasis[0][i]*v[0]+dbasis[1][i]*v[1])*Jdet*gauss->weight;
+                                vx_input->GetInputValue(&vx,gauss);
+                                vy_input->GetInputValue(&vy,gauss);
+                                vel = sqrt(vx*vx+vy*vy);
+                                vx  = vx/(vel+1.e-9);
+                                vy  = vy/(vel+1.e-9);
+                                for(i=0;i<numdof;i++) pe->values[i]+= - weight*(dH[0]*vx+dH[1]*vy)*(dbasis[0][i]*vx+dbasis[1][i]*vy)*Jdet*gauss->weight;
                                 break;
                         case ThicknessAcrossGradientEnum:
                                 weights_input->GetInputValue(&weight, gauss,resp);
-                                vx_input->GetInputValue(&v[0],gauss);
-                                vy_input->GetInputValue(&v[1],gauss);
-                                for(i=0;i<numdof;i++) pe->values[i]+= - weight*(dH[0]*(-v[1])+dH[1]*v[0])*(dbasis[0][i]*(-v[1])+dbasis[1][i]*v[0])*Jdet*gauss->weight;
+                                vx_input->GetInputValue(&vx,gauss);
+                                vy_input->GetInputValue(&vy,gauss);
+                                vel = sqrt(vx*vx+vy*vy);
+                                vx  = vx/(vel+1.e-9);
+                                vy  = vy/(vel+1.e-9);
+                                for(i=0;i<numdof;i++) pe->values[i]+= - weight*(dH[0]*(-vy)+dH[1]*vx)*(dbasis[0][i]*(-vy)+dbasis[1][i]*vx)*Jdet*gauss->weight;
                                 break;
                         default: