Changeset 678

Timestamp:

06/01/09 16:10:31 (16 years ago)

Author:

seroussi

Message:

velocity depth average now working in serial,depth average on the nodes

Location:

issm/trunk/src/c/objects

Files:

• Node.cpp (modified) (1 diff)
• Node.h (modified) (1 diff)
• Penta.cpp (modified) (5 diffs)
• Penta.h (modified) (1 diff)

issm/trunk/src/c/objects/Node.cpp

@@ -491 +491 @@
 }
                 
+#undef __FUNCT__ 
+#define __FUNCT__ "Node::VelocityDepthAverageAtBase"
+void  Node::VelocityDepthAverageAtBase(Vec ug,double* ug_serial){
+
+        /* node data: */
+        int          doflist1;
+        int          dofx,dofy;
+        int          isnodeonsurface;
+        
+        Node* node=NULL;
+        Node* upper_node=NULL;
+        double velocity2[2];
+        double velocity1[2];
+        double velocity_average[2];
+        double sum[2];
+        double z1,z2,dz;
+        double thickness;
+
+        
+        if(onbed==1 & clone==0){
+                        
+                doflist1=GetDofList1();
+
+                /*this node is on the bed. We are going to:
+                 * follow the upper nodes until we reach the surface. At each upper node, we'll grab the 
+                 * velocity for this node, and add it to ug. We'll finish by 
+                 * we grab the velocity. Once we know the velocity, we follow the upper nodes, 
+                 * inserting the same velocity value into ug, until we reach the surface: */
+                sum[0]=0;
+                sum[1]=0;
+                thickness=0;
+
+                /*get dofs for this base node velocity: we know there are two dofs in ug_serial */
+                dofx=2*doflist1;
+                dofy=2*doflist1+1;
+
+                node=this;
+                for(;;){
+
+                        if (node->IsOnSurface())break;
+
+                        doflist1=node->GetDofList1();
+                        
+                        velocity1[0]=ug_serial[2*doflist1];
+                        velocity1[1]=ug_serial[2*doflist1+1];
+                        z1=node->GetZ();
+
+                        upper_node=node->GetUpperNode();
+                        doflist1=upper_node->GetDofList1();
+                
+                        velocity2[0]=ug_serial[2*doflist1];
+                        velocity2[1]=ug_serial[2*doflist1+1];
+                        z2=upper_node->GetZ();
+
+                        dz=(z2-z1);
+                        thickness+=dz;
+                        velocity_average[0]=(velocity1[0]+velocity2[0])/2.0;
+                        velocity_average[1]=(velocity1[1]+velocity2[1])/2.0;
+
+                        sum[0]+=velocity_average[0]*dz;
+                        sum[1]+=velocity_average[1]*dz;
+
+                        /* get next node: */
+                        node=node->GetUpperNode();
+                }
+
+                sum[0]=sum[0]/thickness;
+                sum[1]=sum[1]/thickness;
+
+                /* Plug velocity_average*deltaH/H into base of ug: */
+                VecSetValues(ug,1,&dofx,&sum[0],INSERT_VALUES);
+                VecSetValues(ug,1,&dofy,&sum[1],INSERT_VALUES);
+        }
+}

issm/trunk/src/c/objects/Node.h

@@ -77 +77 @@
                 int   IsOnSurface();
                 void  FreezeDof(int dof);
+                void  VelocityDepthAverageAtBase(Vec ug,double* ug_serial);
 
 };

issm/trunk/src/c/objects/Penta.cpp

@@ -1941 +1941 @@
         /* node data: */
         const int    numgrids=6;
-        const int    numdof=2*numgrids;
-        int          doflist[numdof];
+        int          doflist[numgrids];
+        int          doflist1;
         int          nodedofs[2];
-        int          numberofdofspernode;
         
         Node* node=NULL;
@@ -1953 +1952 @@
         if(onbed==1){
                         
-                GetDofList(&doflist[0],&numberofdofspernode);
+                GetDofList1(&doflist[0]);
 
                 /*this penta is on the bed. For each node on the base of this penta, 
                  * we grab the velocity. Once we know the velocity, we follow the upper nodes, 
-                 * inserting the same velocity value into ug, until we reach the surface: */
+                 * inserting the same velocity value into ug, until we reach the surface: 
+                 * the velocity has two ddls vx and vy*/
                 for(i=0;i<3;i++){
         
@@ -1963 +1963 @@
         
                         /*get velocity for this base node: */
-                        velocity[0]=ug_serial[doflist[numberofdofspernode*i+0]];
-                        velocity[1]=ug_serial[doflist[numberofdofspernode*i+1]];
+                        velocity[0]=ug_serial[2*doflist[i]+0];
+                        velocity[1]=ug_serial[2*doflist[i]+1];
 
                         //go througn all nodes which sit on top of this node, until we reach the surface, 
@@ -1970 +1970 @@
                         for(;;){
 
-                                node->GetDofList(&nodedofs[0],&numberofdofspernode);
+                                doflist1=node->GetDofList1();
+                                nodedofs[0]=2*doflist1;
+                                nodedofs[1]=2*doflist1+1;
                                 VecSetValues(ug,1,&nodedofs[0],&velocity[0],INSERT_VALUES);
                                 VecSetValues(ug,1,&nodedofs[1],&velocity[1],INSERT_VALUES);
@@ -2030 +2032 @@
 
 }
-
-#undef __FUNCT__ 
-#define __FUNCT__ "Penta::VelocityDepthAverageAtBase"
-void  Penta::VelocityDepthAverageAtBase(Vec ug,double* ug_serial){
-
-        /* node data: */
-        const int    numgrids=6;
-        const int    numdof=2*numgrids;
-        int          doflist[numdof];
-        int          nodedofs[2];
-        int          numberofdofspernode;
-        int          dofx,dofy;
-        
-        Node* node=NULL;
-        Node* upper_node=NULL;
-        int   i;
-        double base_velocity[2];
-        double velocity2[2];
-        double velocity1[2];
-        double velocity_average[2];
-        double sum[2];
-        double z1,z2,dz;
-        double thickness;
-                
-        if(onbed==1){
-                        
-                GetDofList(&doflist[0],&numberofdofspernode);
-
-                /*this penta is on the bed. For each node on the base of this penta, we are going to 
-                 * follow the upper nodes until we reach the surface. At each upper node, we'll grab the 
-                 * velocity for this node, and add it to ug. We'll finish by 
-                 * we grab the velocity. Once we know the velocity, we follow the upper nodes, 
-                 * inserting the same velocity value into ug, until we reach the surface: */
-                for(i=0;i<3;i++){
-
-                        //get thickness for this grid
-                        thickness=h[i]; //pick up from the penta element itself.
-        
-                        node=nodes[i]; //base nodes
-        
-                        /*get dofs for this base node velocity: */
-                        dofx=doflist[numberofdofspernode*i+0];
-                        dofy=doflist[numberofdofspernode*i+1];
-
-                        /*first thing, cancel velocity on the base nodes, so that we start with the value 0: */
-                        base_velocity[0]=-ug_serial[doflist[numberofdofspernode*i+0]]; //- to cancel
-                        base_velocity[1]=-ug_serial[doflist[numberofdofspernode*i+1]]; 
-                
-                        VecSetValues(ug,1,&dofx,&base_velocity[0],ADD_VALUES);
-                        VecSetValues(ug,1,&dofy,&base_velocity[1],ADD_VALUES);
-
-                        //go through all nodes which sit on top of this node, until we reach the surface, 
-                        //and plug  deltaV*deltaH/H at base of ug: */
-                        
-                        for(;;){
-
-                                if (node->IsOnSurface())break;
-
-                                node->GetDofList(&nodedofs[0],&numberofdofspernode);
-                                
-                                velocity1[0]=ug_serial[nodedofs[0]];
-                                velocity1[1]=ug_serial[nodedofs[1]];
-                                z1=node->GetZ();
-
-                                upper_node=node->GetUpperNode();
-                                upper_node->GetDofList(&nodedofs[0],&numberofdofspernode);
-                        
-                                velocity2[0]=ug_serial[nodedofs[0]];
-                                velocity2[1]=ug_serial[nodedofs[1]];
-                                z2=upper_node->GetZ();
-
-                                dz=(z2-z1);
-                                velocity_average[0]=(velocity1[0]+velocity2[0])/2.0;
-                                velocity_average[1]=(velocity1[1]+velocity2[1])/2.0;
-
-                                sum[0]=velocity_average[0]*dz/thickness;
-                                sum[1]=velocity_average[1]*dz/thickness;
-
-                                /*Plug velocity_average*deltaH/H into base of ug: */
-                                VecSetValues(ug,1,&dofx,&sum[0],ADD_VALUES);
-                                VecSetValues(ug,1,&dofy,&sum[1],ADD_VALUES);
-
-                                /*get next node: */
-                                node=node->GetUpperNode();
-                        }
-                }
-        }
-}
-
 
 #undef __FUNCT__ 

issm/trunk/src/c/objects/Penta.h

@@ -113 +113 @@
                 void  VelocityExtrude(Vec ug,double* ug_serial);
                 void  ThicknessExtrude(Vec ug,double* ug_serial);
-                void  VelocityDepthAverageAtBase(Vec ug,double* ug_serial);
                 void  SlopeExtrude(Vec sg,double* sg_serial);
                 void  ComputePressure(Vec p_g);