Index: /issm/trunk/src/c/objects/ElementProperties.cpp =================================================================== --- /issm/trunk/src/c/objects/ElementProperties.cpp (revision 3478) +++ /issm/trunk/src/c/objects/ElementProperties.cpp (revision 3479) @@ -34,13 +34,13 @@ this->geothermalflux=NULL; - this->friction_type=0; - this->p=0; - this->q=0; - this->shelf=0; - this->onbed=0; - this->onwater=0; - this->onsurface=0; - this->collapse=0; - this->thermal_steadystate=0; + this->friction_type=UNDEF; + this->p=UNDEF; + this->q=UNDEF; + this->shelf=UNDEF; + this->onbed=UNDEF; + this->onwater=UNDEF; + this->onsurface=UNDEF; + this->collapse=UNDEF; + this->thermal_steadystate=UNDEF; return; } @@ -145,27 +145,27 @@ this->numnodes=elementproperties_numnodes; - if(elementproperties_h)this->h=(double*)xmalloc(this->numnodes*sizeof(double)); - else this->h=NULL; - if(elementproperties_s)this->s=(double*)xmalloc(this->numnodes*sizeof(double)); - else this->s=NULL; - if(elementproperties_b)this->b=(double*)xmalloc(this->numnodes*sizeof(double)); - else this->b=NULL; - if(elementproperties_k)this->k=(double*)xmalloc(this->numnodes*sizeof(double)); - else this->k=NULL; - if(elementproperties_melting)this->melting=(double*)xmalloc(this->numnodes*sizeof(double)); - else this->melting=NULL; - if(elementproperties_accumulation)this->accumulation=(double*)xmalloc(this->numnodes*sizeof(double)); - else this->accumulation=NULL; - if(elementproperties_geothermalflux)this->geothermalflux=(double*)xmalloc(this->numnodes*sizeof(double)); - else this->geothermalflux=NULL; + this->h=(double*)xmalloc(this->numnodes*sizeof(double)); + this->s=(double*)xmalloc(this->numnodes*sizeof(double)); + this->b=(double*)xmalloc(this->numnodes*sizeof(double)); + this->k=(double*)xmalloc(this->numnodes*sizeof(double)); + this->melting=(double*)xmalloc(this->numnodes*sizeof(double)); + this->accumulation=(double*)xmalloc(this->numnodes*sizeof(double)); + this->geothermalflux=(double*)xmalloc(this->numnodes*sizeof(double)); for(i=0;inumnodes;i++){ - if(elementproperties_h)this->h[i]=elementproperties_h[i]; - if(elementproperties_s)this->s[i]=elementproperties_s[i]; - if(elementproperties_b)this->b[i]=elementproperties_b[i]; - if(elementproperties_k)this->k[i]=elementproperties_k[i]; - if(elementproperties_melting)this->melting[i]=elementproperties_melting[i]; - if(elementproperties_accumulation)this->accumulation[i]=elementproperties_accumulation[i]; + if(elementproperties_h) this->h[i]=elementproperties_h[i]; + else this->h[i]=UNDEF; + if(elementproperties_s) this->s[i]=elementproperties_s[i]; + else this->s[i]=UNDEF; + if(elementproperties_b) this->b[i]=elementproperties_b[i]; + else this->b[i]=UNDEF; + if(elementproperties_k) this->k[i]=elementproperties_k[i]; + else this->k[i]=UNDEF; + if(elementproperties_melting) this->melting[i]=elementproperties_melting[i]; + else this->melting[i]=UNDEF; + if(elementproperties_accumulation) this->accumulation[i]=elementproperties_accumulation[i]; + else this->accumulation[i]=UNDEF; if(elementproperties_geothermalflux)this->geothermalflux[i]=elementproperties_geothermalflux[i]; + else this->geothermalflux[i]=UNDEF; } @@ -179,4 +179,53 @@ this->collapse=elementproperties_collapse; this->thermal_steadystate=elementproperties_thermal_steadystate; + + return; +} +/*}}}*/ +/*FUNCTION ElementProperties Init propreties from iomodel{{{1*/ +void ElementProperties::Init(int numberofnodes, int* vertex_ids, int element_id, IoModel* iomodel){ //i is the node index + + int i,k; + + /*Initialize number of nodes*/ + this->numnodes=numberofnodes; + + /*Allocate node properties*/ + if (iomodel->thickness) this->h=(double*)xmalloc(this->numnodes*sizeof(double)); + if (iomodel->surface) this->s=(double*)xmalloc(this->numnodes*sizeof(double)); + if (iomodel->bed) this->b=(double*)xmalloc(this->numnodes*sizeof(double)); + if (iomodel->drag) this->k=(double*)xmalloc(this->numnodes*sizeof(double)); + if (iomodel->melting) this->melting=(double*)xmalloc(this->numnodes*sizeof(double)); + if (iomodel->accumulation) this->accumulation=(double*)xmalloc(this->numnodes*sizeof(double)); + if (iomodel->geothermalflux)this->geothermalflux=(double*)xmalloc(this->numnodes*sizeof(double)); + + /*Get node properties*/ + for(i=0;inumnodes;i++){ + + //get index of the vertex in C indexing + k=vertex_ids[i]-1; + + //Get properties from iomodel + if(iomodel->thickness) this->h[i]=iomodel->thickness[k]; + if(iomodel->surface) this->s[i]=iomodel->surface[k]; + if(iomodel->bed) this->b[i]=iomodel->bed[k]; + if(iomodel->drag) this->k[i]=iomodel->drag[k]; + if(iomodel->melting) this->melting[i]=iomodel->melting[k]; + if(iomodel->accumulation) this->accumulation[i]=iomodel->accumulation[k]; + if(iomodel->geothermalflux)this->geothermalflux[i]=iomodel->geothermalflux[k]; + } + + /*Get element properties from iomodel*/ + if (iomodel->p) this->p=iomodel->p[element_id-1]; + if (iomodel->q) this->q=iomodel->q[element_id-1]; + if (iomodel->elementoniceshelf)this->shelf=(int)iomodel->elementoniceshelf[element_id-1]; + if (iomodel->elementonbed) this->onbed=(int)iomodel->elementonbed[element_id-1]; + if (iomodel->elementonwater) this->onwater=(bool)iomodel->elementonwater[element_id-1]; + if (iomodel->elementonsurface) this->onsurface=(int)iomodel->elementonsurface[element_id-1]; + + /*Get other properties from iomodel*/ + this->friction_type=iomodel->drag_type; + + /*WARNING: collapse and thermal_steadystate are left undefined*/ return; Index: /issm/trunk/src/c/objects/ElementProperties.h =================================================================== --- /issm/trunk/src/c/objects/ElementProperties.h (revision 3478) +++ /issm/trunk/src/c/objects/ElementProperties.h (revision 3479) @@ -5,4 +5,7 @@ #ifndef _ELEMENTPROPERTIES_H_ #define _ELEMENTPROPERTIES_H_ + +/*indefinitions: */ +struct IoModel; class ElementProperties{ @@ -33,7 +36,9 @@ ElementProperties(); ElementProperties(int numnodes, double* h, double* s, double* b, double* k, double* melting, double* accumulation, double* geothermalflux, int friction_type, double p, double q, int shelf, int onbed, bool onwater, int onsurface, int collapse, int thermal_steadystate); - void Init(int numnodes, double* h, double* s, double* b, double* k, double* melting, double* accumulation, double* geothermalflux, int friction_type, double p, double q, int shelf, int onbed, bool onwater, int onsurface, int collapse, int thermal_steadystate); ElementProperties(ElementProperties* properties); ~ElementProperties(); + + void Init(int numnodes, int* vertex_ids, int element_id, IoModel* iomodel); + void Init(int numnodes, double* h, double* s, double* b, double* k, double* melting, double* accumulation, double* geothermalflux, int friction_type, double p, double q, int shelf, int onbed, bool onwater, int onsurface, int collapse, int thermal_steadystate); void Marshall(char** pmarshalled_dataset); Index: /issm/trunk/src/c/objects/Matice.cpp =================================================================== --- /issm/trunk/src/c/objects/Matice.cpp (revision 3478) +++ /issm/trunk/src/c/objects/Matice.cpp (revision 3479) @@ -12,8 +12,8 @@ #include "./Matice.h" #include +#include "../EnumDefinitions/EnumDefinitions.h" #include "../shared/shared.h" #include "../include/macros.h" -#include "../EnumDefinitions/EnumDefinitions.h" - +#include "../include/typedefs.h" /*Object constructors and destructor*/ @@ -51,13 +51,17 @@ matice_mid=i+1; //same as element it is created for - /*B and n: */ - B_avg=0; - for(j=0;jB+((int)*(iomodel->elements+num_vertices*i+j)-1)); - } - B_avg=B_avg/num_vertices; - - matice_B=B_avg; - matice_n=(double)*(iomodel->n+i); + /*Compute B on the element if provided*/ + if (iomodel->B){ + B_avg=0; + for(j=0;jB+((int)*(iomodel->elements+num_vertices*i+j)-1)); + } + B_avg=B_avg/num_vertices; + } + + if (iomodel->B) matice_B=B_avg; + else matice_B=UNDEF; + if (iomodel->n) matice_n=(double)*(iomodel->n+i); + else matice_n=UNDEF; this->Init(matice_mid,matice_B,matice_n); Index: /issm/trunk/src/c/objects/Node.cpp =================================================================== --- /issm/trunk/src/c/objects/Node.cpp (revision 3478) +++ /issm/trunk/src/c/objects/Node.cpp (revision 3479) @@ -803,7 +803,7 @@ /*FUNCTION Node::UpdateFromInputs {{{2*/ void Node::UpdateFromInputs(void* vinputs){ - - ISSMERROR("not used yet!"); - + + /*nothing updated for now*/ + } /*}}}*/ Index: /issm/trunk/src/c/objects/Tria.cpp =================================================================== --- /issm/trunk/src/c/objects/Tria.cpp (revision 3478) +++ /issm/trunk/src/c/objects/Tria.cpp (revision 3479) @@ -67,20 +67,8 @@ int j; - int offset; int tria_node_ids[3]; int tria_matice_id; int tria_matpar_id; int tria_numpar_id; - double tria_h[3]; - double tria_s[3]; - double tria_b[3]; - double tria_k[3]; - double tria_melting[3]; - double tria_accumulation[3]; - int tria_friction_type; - double tria_p,tria_q; - int tria_shelf; - int tria_onwater; - /*id: */ @@ -99,25 +87,5 @@ this->hmatpar.Init(&tria_matpar_id,1); this->hnumpar.Init(&tria_numpar_id,1); - - /*properties: */ - for(j=0;j<3;j++){ - - offset=(int)*(iomodel->elements+3*i+j)-1; //get index of this vertex, in C numbering. - tria_h[j]= *(iomodel->thickness+offset); - tria_s[j]= *(iomodel->surface+offset); - tria_b[j]=*(iomodel->bed+offset); - tria_k[j]=*(iomodel->drag+offset); - tria_melting[j]=*(iomodel->melting+offset); - tria_accumulation[j]=*(iomodel->accumulation+offset); - } - tria_friction_type=(int)iomodel->drag_type; - tria_p=iomodel->p[i]; - tria_q=iomodel->q[i]; - tria_shelf=(int)*(iomodel->elementoniceshelf+i); - tria_onwater=(bool)*(iomodel->elementonwater+i); - - this->properties.Init(3,tria_h, tria_s, tria_b, tria_k, tria_melting, tria_accumulation, NULL, - tria_friction_type, tria_p, tria_q, tria_shelf, UNDEF,tria_onwater, UNDEF,UNDEF,UNDEF); - + this->properties.Init(3,tria_node_ids,this->id,iomodel); } Index: /issm/trunk/src/c/objects/Vertex.cpp =================================================================== --- /issm/trunk/src/c/objects/Vertex.cpp (revision 3478) +++ /issm/trunk/src/c/objects/Vertex.cpp (revision 3479) @@ -188,9 +188,6 @@ Vertex* vertex=NULL; int dof[1]={0}; - double coord[3]; vertex=this; - - coord[0]=this->x; coord[1]=this->y; coord[2]=this->z; /*Recover parameter inputs: */ @@ -198,9 +195,7 @@ /*Update internal data if inputs holds new values: */ - inputs->Recover("x",&coord[0],1,dof,1,(void**)&vertex); - inputs->Recover("y",&coord[1],1,dof,1,(void**)&vertex); - inputs->Recover("z",&coord[2],1,dof,1,(void**)&vertex); - - ISSMERROR("not supported yet!"); + inputs->Recover("x",&x,1,dof,1,(void**)&vertex); + inputs->Recover("y",&y,1,dof,1,(void**)&vertex); + inputs->Recover("z",&z,1,dof,1,(void**)&vertex); } Index: /issm/trunk/src/m/solutions/jpl/control_core.m =================================================================== --- /issm/trunk/src/m/solutions/jpl/control_core.m (revision 3478) +++ /issm/trunk/src/m/solutions/jpl/control_core.m (revision 3479) @@ -45,5 +45,5 @@ %Update inputs in datasets - [model.elements,model.nodes,model.loads,model.materials,model.parameters]=UpdateFromInputs(model.elements,model.nodes,model.loads,model.materials,model.parameters,inputs); + [model.elements,model.nodes,model.vertices,model.loads,model.materials,model.parameters]=UpdateFromInputs(model.elements,model.nodes,model.vertices,model.loads,model.materials,model.parameters,inputs); displaystring(verbose,'\n%s',[' computing gradJ...']); Index: /issm/trunk/src/m/solutions/jpl/diagnostic_core_linear.m =================================================================== --- /issm/trunk/src/m/solutions/jpl/diagnostic_core_linear.m (revision 3478) +++ /issm/trunk/src/m/solutions/jpl/diagnostic_core_linear.m (revision 3479) @@ -9,9 +9,9 @@ %Update inputs in datasets - [m.elements,m.nodes, m.loads,m.materials,m.parameters]=UpdateFromInputs(m.elements,m.nodes, m.loads,m.materials,m.parameters,inputs); + [m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters]=UpdateFromInputs(m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters,inputs); %system matrices - [K_gg, p_g]=SystemMatrices(m.elements,m.nodes,m.loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type); - [K_gg, p_g,kmax]=PenaltySystemMatrices(K_gg,p_g,m.elements,m.nodes,m.loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type); + [K_gg, p_g]=SystemMatrices(m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type); + [K_gg, p_g,kmax]=PenaltySystemMatrices(K_gg,p_g,m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type); %Reduce tangent matrix from g size to f size Index: /issm/trunk/src/m/solutions/jpl/diagnostic_core_nonlinear.m =================================================================== --- /issm/trunk/src/m/solutions/jpl/diagnostic_core_nonlinear.m (revision 3478) +++ /issm/trunk/src/m/solutions/jpl/diagnostic_core_nonlinear.m (revision 3479) @@ -31,11 +31,11 @@ %Update inputs in datasets - [m.elements,m.nodes, loads,m.materials,m.parameters]=UpdateFromInputs(m.elements,m.nodes, loads,m.materials,m.parameters,inputs); + [m.elements,m.nodes,m.vertices, loads,m.materials,m.parameters]=UpdateFromInputs(m.elements,m.nodes,m.vertices, loads,m.materials,m.parameters,inputs); %system matrices - [K_gg_nopenalty , p_g_nopenalty]=SystemMatrices(m.elements,m.nodes,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type); + [K_gg_nopenalty , p_g_nopenalty]=SystemMatrices(m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type); %penalties - [K_gg , p_g, kmax]=PenaltySystemMatrices(K_gg_nopenalty,p_g_nopenalty,m.elements,m.nodes,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type); + [K_gg , p_g, kmax]=PenaltySystemMatrices(K_gg_nopenalty,p_g_nopenalty,m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type); %Reduce tangent matrix from g size to f size @@ -84,5 +84,5 @@ inputs=add(inputs,'velocity',u_g,'doublevec',m.parameters.numberofdofspernode,m.parameters.numberofnodes); m.parameters.kflag=1; m.parameters.pflag=0; - [K_gg, p_g]=SystemMatrices(m.elements,m.nodes,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type); + [K_gg, p_g]=SystemMatrices(m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type); [K_ff, K_fs] = Reducematrixfromgtof( K_gg, m.Gmn, m.nodesets); varargout(1)={K_ff}; Index: /issm/trunk/src/m/solutions/jpl/thermal_core_nonlinear.m =================================================================== --- /issm/trunk/src/m/solutions/jpl/thermal_core_nonlinear.m (revision 3478) +++ /issm/trunk/src/m/solutions/jpl/thermal_core_nonlinear.m (revision 3479) @@ -23,5 +23,5 @@ %Update inputs in datasets - [m.elements,m.nodes, loads,m.materials,m.parameters]=UpdateFromInputs(m.elements,m.nodes, loads,m.materials,m.parameters,inputs); + [m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters]=UpdateFromInputs(m.elements,m.nodes,m.vertices, loads,m.materials,m.parameters,inputs); %system matrices @@ -29,13 +29,13 @@ if count==1 displaystring(m.parameters.verbose,'%s',[' system matrices']); - [K_gg_nopenalty, p_g_nopenalty]=SystemMatrices(m.elements,m.nodes,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type); + [K_gg_nopenalty, p_g_nopenalty]=SystemMatrices(m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type); end displaystring(m.parameters.verbose,'%s',[' penalty system matrices']); - [K_gg , p_g, melting_offset]=PenaltySystemMatrices(K_gg_nopenalty,p_g_nopenalty,m.elements,m.nodes,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type); + [K_gg , p_g, melting_offset]=PenaltySystemMatrices(K_gg_nopenalty,p_g_nopenalty,m.elements,m.nodes,m.verticesloads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type); else displaystring(m.parameters.verbose,'%s',[' system matrices']); - [K_gg , p_g]=SystemMatrices(m.elements,m.nodes,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type); + [K_gg , p_g]=SystemMatrices(m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type); displaystring(m.parameters.verbose,'%s',[' penalty system matrices']); - [K_gg , p_g, melting_offset]=PenaltySystemMatrices(K_gg,p_g,m.elements,m.nodes,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type); + [K_gg , p_g, melting_offset]=PenaltySystemMatrices(K_gg,p_g,m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type); end @@ -57,5 +57,5 @@ inputs=add(inputs,'temperature',t_g,'doublevec',m.parameters.numberofdofspernode,m.parameters.numberofnodes); displaystring(m.parameters.verbose,'%s',[' update inputs']); - [m.elements,m.nodes, loads,m.materials,m.parameters]=UpdateFromInputs(m.elements,m.nodes, loads,m.materials,m.parameters,inputs); + [m.elements,m.nodes,m.vertices, loads,m.materials,m.parameters]=UpdateFromInputs(m.elements,m.nodes,m.vertices, loads,m.materials,m.parameters,inputs); %penalty constraints Index: /issm/trunk/src/mex/PenaltySystemMatrices/PenaltySystemMatrices.cpp =================================================================== --- /issm/trunk/src/mex/PenaltySystemMatrices/PenaltySystemMatrices.cpp (revision 3478) +++ /issm/trunk/src/mex/PenaltySystemMatrices/PenaltySystemMatrices.cpp (revision 3479) @@ -79,5 +79,5 @@ { _printf_("\n"); - _printf_(" usage: [Kgg,pg] = %s(Kggin,pgin,eleemnts,nodes,loads,materials,params,inputs,analysis_type,sub_analysis_type);\n",__FUNCT__); + _printf_(" usage: [Kgg,pg] = %s(Kggin,pgin,elements,nodes,vertices,loads,materials,params,inputs,analysis_type,sub_analysis_type);\n",__FUNCT__); _printf_("\n"); } Index: /issm/trunk/src/mex/SystemMatrices/SystemMatrices.cpp =================================================================== --- /issm/trunk/src/mex/SystemMatrices/SystemMatrices.cpp (revision 3478) +++ /issm/trunk/src/mex/SystemMatrices/SystemMatrices.cpp (revision 3479) @@ -81,5 +81,5 @@ { _printf_("\n"); - _printf_(" usage: [Kgg,pg] = %s(eleemnts,nodes,loads,materials,params,inputs,analysis_type);\n",__FUNCT__); + _printf_(" usage: [Kgg,pg] = %s(elements,nodes,vertices,loads,materials,params,inputs,analysis_type);\n",__FUNCT__); _printf_("\n"); } Index: /issm/trunk/src/mex/UpdateFromInputs/UpdateFromInputs.cpp =================================================================== --- /issm/trunk/src/mex/UpdateFromInputs/UpdateFromInputs.cpp (revision 3478) +++ /issm/trunk/src/mex/UpdateFromInputs/UpdateFromInputs.cpp (revision 3479) @@ -64,5 +64,5 @@ { _printf_("\n"); - _printf_(" usage: [elements,nodes,loads, materials, parameters] = %s(elements,nodes,loads, materials,parameters,inputs);\n",__FUNCT__); + _printf_(" usage: [elements,nodes,vertices,loads,materials,parameters] = %s(elements,nodes,vertices,loads,materials,parameters,inputs);\n",__FUNCT__); _printf_("\n"); }