Index: /issm/trunk/src/c/EnumDefinitions/EnumDefinitions.cpp =================================================================== --- /issm/trunk/src/c/EnumDefinitions/EnumDefinitions.cpp (revision 307) +++ /issm/trunk/src/c/EnumDefinitions/EnumDefinitions.cpp (revision 308) @@ -42,4 +42,6 @@ int TriaEnum(void){ return 101; } int PentaEnum(void){ return 102; } +int SingEnum(void){ return 103; } +int BeamEnum(void){ return 104; } /*Loads: */ Index: /issm/trunk/src/c/EnumDefinitions/EnumDefinitions.h =================================================================== --- /issm/trunk/src/c/EnumDefinitions/EnumDefinitions.h (revision 307) +++ /issm/trunk/src/c/EnumDefinitions/EnumDefinitions.h (revision 308) @@ -10,4 +10,6 @@ int NodeEnum(void); int TriaEnum(void); +int SingEnum(void); +int BeamEnum(void); int PentaEnum(void); int MaterialEnum(void); Index: /issm/trunk/src/c/Makefile.am =================================================================== --- /issm/trunk/src/c/Makefile.am (revision 307) +++ /issm/trunk/src/c/Makefile.am (revision 308) @@ -36,4 +36,8 @@ ./objects/Tria.h\ ./objects/Tria.cpp\ + ./objects/Sing.h\ + ./objects/Sing.cpp\ + ./objects/Beam.h\ + ./objects/Beam.cpp\ ./objects/Penta.h\ ./objects/Penta.cpp\ @@ -273,4 +277,8 @@ ./objects/Tria.h\ ./objects/Tria.cpp\ + ./objects/Sing.h\ + ./objects/Sing.cpp\ + ./objects/Beam.h\ + ./objects/Beam.cpp\ ./objects/Penta.h\ ./objects/Penta.cpp\ Index: /issm/trunk/src/c/ModelProcessorx/DiagnosticHutter/CreateConstraintsDiagnosticHutter.cpp =================================================================== --- /issm/trunk/src/c/ModelProcessorx/DiagnosticHutter/CreateConstraintsDiagnosticHutter.cpp (revision 307) +++ /issm/trunk/src/c/ModelProcessorx/DiagnosticHutter/CreateConstraintsDiagnosticHutter.cpp (revision 308) @@ -16,5 +16,5 @@ void CreateConstraintsDiagnosticHutter(DataSet** pconstraints, Model* model,ConstDataHandle model_handle){ - int i; + int i,j; int count; @@ -22,4 +22,5 @@ Spc* spc = NULL; + Rgb* rgb = NULL; /*spc intermediary data: */ @@ -28,5 +29,11 @@ int spc_dof; double spc_value; - + + /*rgb constructor data: */ + int rgb_id; + int rgb_dof; + int rgb_nodeid1; + int rgb_nodeid2; + /*Create constraints: */ constraints = new DataSet(ConstraintsEnum()); @@ -90,6 +97,6 @@ rgb_id=count; rgb_dof=1; - rgb_nodeid1=*(model->penalties+(model->numlayers-1)*i+j); - rgb_nodeid1=*(model->penalties+(model->numlayers-1)*i+j+1); + rgb_nodeid1=(int)*(model->penalties+(model->numlayers-1)*i+j); + rgb_nodeid1=(int)*(model->penalties+(model->numlayers-1)*i+j+1); rgb = new Rgb(rgb_id,rgb_nodeid1,rgb_nodeid2,rgb_dof); @@ -100,6 +107,6 @@ rgb_id=count; rgb_dof=2; - rgb_nodeid1=*(model->penalties+(model->numlayers-1)*i+j); - rgb_nodeid1=*(model->penalties+(model->numlayers-1)*i+j+1); + rgb_nodeid1=(int)*(model->penalties+(model->numlayers-1)*i+j); + rgb_nodeid1=(int)*(model->penalties+(model->numlayers-1)*i+j+1); rgb = new Rgb(rgb_id,rgb_nodeid1,rgb_nodeid2,rgb_dof); Index: /issm/trunk/src/c/ModelProcessorx/DiagnosticHutter/CreateElementsNodesAndMaterialsDiagnosticHutter.cpp =================================================================== --- /issm/trunk/src/c/ModelProcessorx/DiagnosticHutter/CreateElementsNodesAndMaterialsDiagnosticHutter.cpp (revision 307) +++ /issm/trunk/src/c/ModelProcessorx/DiagnosticHutter/CreateElementsNodesAndMaterialsDiagnosticHutter.cpp (revision 308) @@ -17,4 +17,8 @@ void CreateElementsNodesAndMaterialsDiagnosticHutter(DataSet** pelements,DataSet** pnodes, DataSet** pmaterials, Model* model,ConstDataHandle model_handle){ + + int i,k; + extern int my_rank; + extern int num_procs; /*DataSets: */ @@ -68,4 +72,5 @@ double beam_b[2]; double beam_k[2]; + int beam_onbed; /*matpar constructor input: */ @@ -139,5 +144,6 @@ ModelFetchData((void**)&model->surface,NULL,NULL,model_handle,"surface","Matrix","Mat"); ModelFetchData((void**)&model->gridonsurface,NULL,NULL,model_handle,"gridonsurface","Matrix","Mat"); - ModelFetchData((void**)&model->uppergrids,NULL,NULL,model_handle,"uppergrids","Matrix","Mat"); + ModelFetchData((void**)&model->gridonbed,NULL,NULL,model_handle,"gridonbed","Matrix","Mat"); + ModelFetchData((void**)&model->uppernodes,NULL,NULL,model_handle,"uppergrids","Matrix","Mat"); ModelFetchData((void**)&model->drag,NULL,NULL,model_handle,"drag","Matrix","Mat"); ModelFetchData((void**)&model->B,NULL,NULL,model_handle,"B","Matrix","Mat"); @@ -204,16 +210,18 @@ beam_mparid=model->numberofnodes+1;//refers to the corresponding matpar property card beam_g[0]=i+1; - beam_g[1]=model->uppergrids[i]; //grid that lays right on top + beam_g[1]=(int)model->uppernodes[i]; //grid that lays right on top beam_h[0]=model->thickness[i]; - beam_h[1]=model->thickness[model->uppergrids[i]-1]; + beam_h[1]=model->thickness[(int)(model->uppernodes[i]-1)]; beam_s[0]=model->surface[i]; - beam_s[1]=model->surface[model->uppergrids[i]-1]; + beam_s[1]=model->surface[(int)(model->uppernodes[i]-1)]; beam_b[0]=model->bed[i]; - beam_b[1]=model->bed[model->uppergrids[i]-1]; + beam_b[1]=model->bed[(int)(model->uppernodes[i]-1)]; beam_k[0]=model->drag[i]; - beam_k[1]=model->drag[model->uppergrids[i]-1]; + beam_k[1]=model->drag[(int)(model->uppernodes[i]-1)]; + + if (model->gridonbed[i] ) beam_onbed=1; /*Create beam element ubeam its constructor:*/ - beam=new Beam(beam_id, beam_mid, beam_mparid, beam_g, beam_h, beam_s,beam_b,beam_k); + beam=new Beam(beam_id, beam_mid, beam_mparid, beam_g, beam_h, beam_s,beam_b,beam_k,beam_onbed); /*Add tria element to elements dataset: */ @@ -247,5 +255,5 @@ xfree((void**)&model->surface); xfree((void**)&model->gridonsurface); - xfree((void**)&model->uppergrids); + xfree((void**)&model->uppernodes); xfree((void**)&model->drag); xfree((void**)&model->B); Index: /issm/trunk/src/c/ModelProcessorx/DiagnosticHutter/CreateLoadsDiagnosticHutter.cpp =================================================================== --- /issm/trunk/src/c/ModelProcessorx/DiagnosticHutter/CreateLoadsDiagnosticHutter.cpp (revision 307) +++ /issm/trunk/src/c/ModelProcessorx/DiagnosticHutter/CreateLoadsDiagnosticHutter.cpp (revision 308) @@ -30,4 +30,2 @@ } - - Index: /issm/trunk/src/c/ModelProcessorx/Model.cpp =================================================================== --- /issm/trunk/src/c/ModelProcessorx/Model.cpp (revision 307) +++ /issm/trunk/src/c/ModelProcessorx/Model.cpp (revision 308) @@ -49,4 +49,5 @@ model->numlayers=0; model->uppernodes=NULL; + model->gridonhutter=NULL; model->vx_obs=NULL; @@ -183,4 +184,5 @@ xfree((void**)&model->elements); xfree((void**)&model->elements_type); + xfree((void**)&model->gridonhutter); if (strcmp(model->meshtype,"3d")==0){ xfree((void**)&model->elements2d); Index: /issm/trunk/src/c/ModelProcessorx/Model.h =================================================================== --- /issm/trunk/src/c/ModelProcessorx/Model.h (revision 307) +++ /issm/trunk/src/c/ModelProcessorx/Model.h (revision 308) @@ -39,4 +39,5 @@ int ismacayealpattyn; int isstokes; + double* gridonhutter; /*results: */ Index: /issm/trunk/src/c/objects/Beam.cpp =================================================================== --- /issm/trunk/src/c/objects/Beam.cpp (revision 308) +++ /issm/trunk/src/c/objects/Beam.cpp (revision 308) @@ -0,0 +1,608 @@ +/*!\file Beam.c + * \brief: implementation of the Beam object + */ + +#ifdef HAVE_CONFIG_H + #include "config.h" +#else +#error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!" +#endif + +#include "stdio.h" +#include "./Beam.h" +#include +#include "../EnumDefinitions/EnumDefinitions.h" +#include "../shared/shared.h" +#include "../DataSet/DataSet.h" +#include "../include/typedefs.h" + + +Beam::Beam(){ + return; +} + +Beam::~Beam(){ + return; +} + +Beam::Beam(int beam_id, int beam_mid, int beam_mparid, int beam_g[2], double beam_h[2], double beam_s[2],double beam_b[2],double beam_k[2],int beam_onbed){ + + int i; + + id=beam_id; + mid=beam_mid; + mparid=beam_mparid; + for(i=0;i<2;i++){ + node_ids[i]=beam_g[i]; + node_offsets[i]=UNDEF; + nodes[i]=NULL; + h[i]=beam_h[i]; + s[i]=beam_s[i]; + b[i]=beam_b[i]; + k[i]=beam_k[i]; + } + matice=NULL; + matice_offset=UNDEF; + matpar=NULL; + matpar_offset=UNDEF; + onbed=beam_onbed; + + return; +} + +#undef __FUNCT__ +#define __FUNCT__ "Beam::Echo" + +void Beam::Echo(void){ + + printf("Beam:\n"); + printf(" id: %i\n",id); + printf(" mid: %i\n",mid); + printf(" mparid: %i\n",mparid); + printf(" node_ids=[%i,%i]\n",node_ids[0],node_ids[1]); + printf(" node_offsets=[%i,%i]\n",node_offsets[0],node_offsets[1]); + printf(" matice_offset=%i\n",matice_offset); + printf(" matpar_offset=%i\n",matpar_offset); + printf(" h=[%g,%g]\n",h[0],h[1]); + printf(" s=[%g,%g]\n",s[0],s[1]); + printf(" b=[%g,%g]\n",b[0],b[1]); + printf(" k=[%g,%g]\n",k[0],k[1]); + printf(" onbed=%i\n",onbed); + printf(" nodes: \n"); + if(nodes[0])nodes[0]->Echo(); + if(nodes[1])nodes[1]->Echo(); + if(matice)matice->Echo(); + if(matpar)matpar->Echo(); + + return; +} + +void Beam::Marshall(char** pmarshalled_dataset){ + + char* marshalled_dataset=NULL; + int enum_type=0; + + /*recover marshalled_dataset: */ + marshalled_dataset=*pmarshalled_dataset; + + /*get enum type of Beam: */ + enum_type=BeamEnum(); + + /*marshall enum: */ + memcpy(marshalled_dataset,&enum_type,sizeof(enum_type));marshalled_dataset+=sizeof(enum_type); + + /*marshall Beam data: */ + memcpy(marshalled_dataset,&id,sizeof(id));marshalled_dataset+=sizeof(id); + memcpy(marshalled_dataset,&mid,sizeof(mid));marshalled_dataset+=sizeof(mid); + memcpy(marshalled_dataset,&mparid,sizeof(mparid));marshalled_dataset+=sizeof(mparid); + memcpy(marshalled_dataset,&node_ids,sizeof(node_ids));marshalled_dataset+=sizeof(node_ids); + memcpy(marshalled_dataset,&nodes,sizeof(nodes));marshalled_dataset+=sizeof(nodes); + memcpy(marshalled_dataset,&node_offsets,sizeof(node_offsets));marshalled_dataset+=sizeof(node_offsets); + memcpy(marshalled_dataset,&matice,sizeof(matice));marshalled_dataset+=sizeof(matice); + memcpy(marshalled_dataset,&matice_offset,sizeof(matice_offset));marshalled_dataset+=sizeof(matice_offset); + memcpy(marshalled_dataset,&matpar,sizeof(matpar));marshalled_dataset+=sizeof(matpar); + memcpy(marshalled_dataset,&matpar_offset,sizeof(matpar_offset));marshalled_dataset+=sizeof(matpar_offset); + memcpy(marshalled_dataset,&h,sizeof(h));marshalled_dataset+=sizeof(h); + memcpy(marshalled_dataset,&s,sizeof(s));marshalled_dataset+=sizeof(s); + memcpy(marshalled_dataset,&b,sizeof(b));marshalled_dataset+=sizeof(b); + memcpy(marshalled_dataset,&k,sizeof(k));marshalled_dataset+=sizeof(k); + memcpy(marshalled_dataset,&onbed,sizeof(onbed));marshalled_dataset+=sizeof(onbed); + + *pmarshalled_dataset=marshalled_dataset; + return; +} + +int Beam::MarshallSize(){ + return sizeof(id) + +sizeof(mid) + +sizeof(mparid) + +sizeof(node_ids) + +sizeof(nodes) + +sizeof(node_offsets) + +sizeof(matice) + +sizeof(matice_offset) + +sizeof(matpar) + +sizeof(matpar_offset) + +sizeof(h) + +sizeof(s) + +sizeof(b) + +sizeof(k) + +sizeof(onbed) + +sizeof(int); //sizeof(int) for enum type +} + +char* Beam::GetName(void){ + return "beam"; +} + +void Beam::Demarshall(char** pmarshalled_dataset){ + + char* marshalled_dataset=NULL; + int i; + + /*recover marshalled_dataset: */ + marshalled_dataset=*pmarshalled_dataset; + + /*this time, no need to get enum type, the pointer directly points to the beginning of the + *object data (thanks to DataSet::Demarshall):*/ + + memcpy(&id,marshalled_dataset,sizeof(id));marshalled_dataset+=sizeof(id); + memcpy(&mid,marshalled_dataset,sizeof(mid));marshalled_dataset+=sizeof(mid); + memcpy(&mparid,marshalled_dataset,sizeof(mparid));marshalled_dataset+=sizeof(mparid); + memcpy(&node_ids,marshalled_dataset,sizeof(node_ids));marshalled_dataset+=sizeof(node_ids); + memcpy(&nodes,marshalled_dataset,sizeof(nodes));marshalled_dataset+=sizeof(nodes); + memcpy(&node_offsets,marshalled_dataset,sizeof(node_offsets));marshalled_dataset+=sizeof(node_offsets); + memcpy(&matice,marshalled_dataset,sizeof(matice));marshalled_dataset+=sizeof(matice); + memcpy(&matice_offset,marshalled_dataset,sizeof(matice_offset));marshalled_dataset+=sizeof(matice_offset); + memcpy(&matpar,marshalled_dataset,sizeof(matpar));marshalled_dataset+=sizeof(matpar); + memcpy(&matpar_offset,marshalled_dataset,sizeof(matpar_offset));marshalled_dataset+=sizeof(matpar_offset); + memcpy(&h,marshalled_dataset,sizeof(h));marshalled_dataset+=sizeof(h); + memcpy(&s,marshalled_dataset,sizeof(s));marshalled_dataset+=sizeof(s); + memcpy(&b,marshalled_dataset,sizeof(b));marshalled_dataset+=sizeof(b); + memcpy(&k,marshalled_dataset,sizeof(k));marshalled_dataset+=sizeof(k); + memcpy(&onbed,marshalled_dataset,sizeof(onbed));marshalled_dataset+=sizeof(onbed); + + /*nodes and materials are not pointing to correct objects anymore:*/ + for(i=0;i<2;i++)nodes[i]=NULL; + matice=NULL; + matpar=NULL; + + /*return: */ + *pmarshalled_dataset=marshalled_dataset; + return; +} +int Beam::Enum(void){ + + return BeamEnum(); + +} +int Beam::GetId(void){ return id; } + +int Beam::MyRank(void){ + extern int my_rank; + return my_rank; +} + + +#undef __FUNCT__ +#define __FUNCT__ "Beam::Configure" +void Beam::Configure(void* ploadsin,void* pnodesin,void* pmaterialsin){ + + int i; + + DataSet* loadsin=NULL; + DataSet* nodesin=NULL; + DataSet* materialsin=NULL; + + /*Recover pointers :*/ + loadsin=(DataSet*)ploadsin; + nodesin=(DataSet*)pnodesin; + materialsin=(DataSet*)pmaterialsin; + + /*Link this element with its nodes, ie find pointers to the nodes in the nodes dataset.: */ + ResolvePointers((Object**)nodes,node_ids,node_offsets,2,nodesin); + + /*Same for materials: */ + ResolvePointers((Object**)&matice,&mid,&matice_offset,1,materialsin); + ResolvePointers((Object**)&matpar,&mparid,&matpar_offset,1,materialsin); + +} + +#undef __FUNCT__ +#define __FUNCT__ "Beam::CreateKMatrix" + +void Beam::CreateKMatrix(Mat Kgg,void* inputs,int analysis_type){ + + /*Just branch to the correct element stiffness matrix generator, according to the type of analysis we are carrying out: */ + if (analysis_type==DiagnosticHutterAnalysisEnum()) { + + CreateKMatrixDiagnosticHutter( Kgg,inputs,analysis_type); + + } + else{ + throw ErrorException(__FUNCT__,exprintf("%s%i%s\n","analysis: ",analysis_type," not supported yet")); + } + +} + + +#undef __FUNCT__ +#define __FUNCT__ "Beam::CreateKMatrixDiagnosticHutter" + +void Beam::CreateKMatrixDiagnosticHutter(Mat Kgg,void* vinputs,int analysis_type){ + + + const int numgrids=2; + const int NDOF2=2; + const int numdofs=NDOF2*numgrids; + int doflist[numdofs]; + double Ke_gg[numdofs][numdofs]={{0,0,0,0},{0,0,0,0},{0,0,0,0},{0,0,0,0}}; + int numberofdofspernode; + + GetDofList(&doflist[0],&numberofdofspernode); + + if (onbed){ + Ke_gg[0][0]=1; + Ke_gg[1][1]=1; + Ke_gg[2][0]=-1; + Ke_gg[2][2]=1; + Ke_gg[3][1]=-1; + Ke_gg[3][3]=1; + } + else{ + Ke_gg[2][0]=-1; + Ke_gg[2][2]=1; + Ke_gg[3][1]=-1; + Ke_gg[3][3]=1; + } + + /*Add Ke_gg to global matrix Kgg: */ + MatSetValues(Kgg,numdofs,doflist,numdofs,doflist,(const double*)Ke_gg,ADD_VALUES); + +} + + +#undef __FUNCT__ +#define __FUNCT__ "Beam::CreatePVector" +void Beam::CreatePVector(Vec pg,void* inputs,int analysis_type){ + + /*Just branch to the correct load generator, according to the type of analysis we are carrying out: */ + if (analysis_type==DiagnosticHutterAnalysisEnum()) { + CreatePVectorDiagnosticHutter( pg,inputs,analysis_type); + } + else{ + throw ErrorException(__FUNCT__,exprintf("%s%i%s"," analysis ",analysis_type," not supported yet")); + } + +} + + +#undef __FUNCT__ +#define __FUNCT__ "Beam::CreatePVectorDiagnosticHutter" + +void Beam::CreatePVectorDiagnosticHutter( Vec pg, void* vinputs, int analysis_type){ + + int i,j,k; + + const int numgrids=2; + const int NDOF2=2; + const int numdofs=NDOF2*numgrids; + int doflist[numdofs]; + double pe_g[4]={0,0,0,0}; + double pe_g_gaussian[4]={0,0,0,0}; + int found=0; + int dofs[1]={0}; + double xyz_list[numgrids][3]; + double z_list[numgrids]; + double constant_part; + int numberofdofspernode; + + /*gaussian points: */ + int num_gauss; + double* segment_gauss_coord=NULL; + double gauss_coord; + double* gauss_weights=NULL; + double gauss_weight; + int ig; + double l1l2[2]; + double slope[2]; + double slope2; + + double z_g; + double rho_ice,gravity,n,B; + double Jdet; + double ub,vb; + + ParameterInputs* inputs=NULL; + + /*recover pointers: */ + inputs=(ParameterInputs*)vinputs; + + /*recover doflist: */ + GetDofList(&doflist[0],&numberofdofspernode); + + /*recover parameters: */ + rho_ice=matpar->GetRhoIce(); + gravity=matpar->GetG(); + n=matice->GetN(); + B=matice->GetB(); + + //recover extra inputs + found=inputs->Recover("surfaceslopex",&slope[0],1,dofs,1,(void**)&nodes[0]); //only recover from fist node, second node will have the same slope + if(!found)throw ErrorException(__FUNCT__," surfaceslopex missing from inputs!"); + + found=inputs->Recover("surfaceslopey",&slope[1],1,dofs,1,(void**)&nodes[0]);//only recover from fist node, second node will have the same slope + if(!found)throw ErrorException(__FUNCT__," surfaceslopey missing from inputs!"); + + //Get all element grid data: + GetElementNodeData( &xyz_list[0][0], nodes, numgrids); + for(i=0;iRecover("thickness",&h[0],1,dofs,2,(void**)nodes); + inputs->Recover("surface",&s[0],1,dofs,2,(void**)nodes); + inputs->Recover("bed",&b[0],1,dofs,2,(void**)nodes); + inputs->Recover("drag",&k[0],1,dofs,2,(void**)nodes); + + //Update material if necessary + if(inputs->Recover("temperature_average",&temperature_list[0],1,dofs,2,(void**)nodes)){ + temperature_average=(temperature_list[0]+temperature_list[1])/2.0; + B_average=Paterson(temperature_average); + matice->SetB(B_average); + } + if(inputs->Recover("temperature",&temperature_list[0],1,dofs,2,(void**)nodes)){ + temperature_average=(temperature_list[0]+temperature_list[1])/2.0; + B_average=Paterson(temperature_average); + matice->SetB(B_average); + } + + if(inputs->Recover("B",&B_list[0],1,dofs,2,(void**)nodes)){ + B_average=(B_list[0]+B_list[1])/2.0; + matice->SetB(B_average); + } + +} + +void Beam::GetDofList(int* doflist,int* pnumberofdofspernode){ + + int i,j; + int doflist_per_node[MAXDOFSPERNODE]; + int numberofdofspernode; + + for(i=0;i<2;i++){ + nodes[i]->GetDofList(&doflist_per_node[0],&numberofdofspernode); + for(j=0;jGetDofList1(); + } + +} + + +void* Beam::GetMatPar(){ + return matpar; +} + + +void Beam::GetNodes(void** vpnodes){ + int i; + Node** pnodes=(Node**)vpnodes; + + for(i=0;i<3;i++){ + pnodes[i]=nodes[i]; + } +} + + +Object* Beam::copy() { + + return new Beam(*this); + +} + + +#undef __FUNCT__ +#define __FUNCT__ "Beam::NodeConfiguration" +void Beam::NodeConfiguration(int* beam_node_ids,Node* beam_nodes[2],int* beam_node_offsets){ + + int i; + for(i=0;i<2;i++){ + node_ids[i]=beam_node_ids[i]; + nodes[i]=beam_nodes[i]; + node_offsets[i]=beam_node_offsets[i]; + } + +} +#undef __FUNCT__ +#define __FUNCT__ "Beam::MaticeConfiguration" +void Beam::MaticeConfiguration(Matice* beam_matice,int beam_matice_offset){ + matice=beam_matice; + matice_offset=beam_matice_offset; +} + +#undef __FUNCT__ +#define __FUNCT__ "Beam::MatparConfiguration" +void Beam::MatparConfiguration(Matpar* beam_matpar,int beam_matpar_offset){ + + matpar=beam_matpar; + matpar_offset=beam_matpar_offset; + +} + + + +#undef __FUNCT__ +#define __FUNCT__ "Beam::GetShelf" +int Beam::GetShelf(){ + throw ErrorException(__FUNCT__," not supported yet!"); +} +#undef __FUNCT__ +#define __FUNCT__ "Beam::GetOnBed" +int Beam::GetOnBed(){ + throw ErrorException(__FUNCT__," not supported yet!"); +} + + +#undef __FUNCT__ +#define __FUNCT__ "Beam::GetBedList" +void Beam::GetBedList(double*){ + throw ErrorException(__FUNCT__," not supported yet!"); +} + +#undef __FUNCT__ +#define __FUNCT__ "Beam::Du" +void Beam::Du(_p_Vec*, double*, double*, void*, int){ + throw ErrorException(__FUNCT__," not supported yet!"); +} +#undef __FUNCT__ +#define __FUNCT__ "Beam::Gradj" +void Beam::Gradj(_p_Vec*, double*, double*, void*, int, char*){ + throw ErrorException(__FUNCT__," not supported yet!"); +} +#undef __FUNCT__ +#define __FUNCT__ "Beam::GradjDrag" +void Beam::GradjDrag(_p_Vec*, double*, double*, void*, int){ + throw ErrorException(__FUNCT__," not supported yet!"); +} +#undef __FUNCT__ +#define __FUNCT__ "Beam::GradjB" +void Beam::GradjB(_p_Vec*, double*, double*, void*, int){ + throw ErrorException(__FUNCT__," not supported yet!"); +} +#undef __FUNCT__ +#define __FUNCT__ "Beam::Misfit" +double Beam::Misfit(double*, double*, void*, int){ + throw ErrorException(__FUNCT__," not supported yet!"); +} + +#undef __FUNCT__ +#define __FUNCT__ "Beam::GetThicknessList" +void Beam::GetThicknessList(double* thickness_list){ + throw ErrorException(__FUNCT__," not supported yet!"); +} + +#undef __FUNCT__ +#define __FUNCT__ "Beam::ComputePressure" +void Beam::ComputePressure(Vec p_g){ + throw ErrorException(__FUNCT__," not supported yet!"); +} + +#undef __FUNCT__ +#define __FUNCT__ "Beam::GetParameterValue" +void Beam::GetParameterValue(double* pvalue, double* value_list,double gauss_coord){ + + double l1l2[2]; + + GetNodalFunctions(&l1l2[0],gauss_coord); + + *pvalue=l1l2[0]*value_list[0]+l1l2[1]*value_list[1]; +} + + +#undef __FUNCT__ +#define __FUNCT__ "Beam::GetNodalFunctions" +void Beam::GetNodalFunctions(double* l1l2, double gauss_coord){ + + /*This routine returns the values of the nodal functions at the gaussian point.*/ + + /*First nodal function: */ + l1l2[0]=.5*gauss_coord+.5; + + /*Second nodal function: */ + l1l2[1]=-.5*gauss_coord+.5; +} + +#undef __FUNCT__ +#define __FUNCT__ "Beam::GetJacobianDeterminant" +void Beam::GetJacobianDeterminant(double* pJdet,double* z_list, double gauss_coord){ + + + double Jdet; + + Jdet=1.0/2.0*(z_list[1]-z_list[0]); + + if(Jdet<0){ + throw ErrorException(__FUNCT__," negative jacobian determinant!"); + } + + *pJdet=Jdet; +} + Index: /issm/trunk/src/c/objects/Beam.h =================================================================== --- /issm/trunk/src/c/objects/Beam.h (revision 308) +++ /issm/trunk/src/c/objects/Beam.h (revision 308) @@ -0,0 +1,90 @@ +/*! \file Beam.h + * \brief: header file for beam object + */ + +#ifndef _BEAM_H_ +#define _BEAM_H_ + + +#include "./Element.h" +#include "./Node.h" +#include "./Matice.h" +#include "./Matpar.h" +#include "./ParameterInputs.h" + +class Beam: public Element{ + + private: + + /*ids:*/ + int id; + + /*nodes: */ + int node_ids[2]; //node ids + Node* nodes[2]; //node pointers + int node_offsets[2]; //node offsets in nodes dataset + + /*materials: */ + int mid; + Matice* matice; + int matice_offset; + + int mparid; + Matpar* matpar; + int matpar_offset; + + double h[2]; + double s[2]; + double b[2]; + double k[2]; + + int onbed; + + public: + + Beam(); + Beam(int beam_id, int beam_mid, int beam_mparid, int beam_g[2], double beam_h[2], double beam_s[2],double beam_b[2],double beam_k[2],int beam_onbed); + ~Beam(); + + void Echo(); + void Marshall(char** pmarshalled_dataset); + int MarshallSize(); + char* GetName(); + void Demarshall(char** pmarshalled_dataset); + int Enum(); + int GetId(); + int MyRank(); + void Configure(void* loads,void* nodes,void* materials); + void CreateKMatrix(Mat Kgg,void* inputs,int analysis_type); + void CreatePVector(Vec pg, void* inputs, int analysis_type); + void UpdateFromInputs(void* inputs); + void GetDofList(int* doflist,int* pnumberofdofs); + void GetDofList1(int* doflist); + void CreateKMatrixDiagnosticHutter(Mat Kgg,void* inputs,int analysis_type); + void GetParameterValue(double* pp, double* plist, double* gauss_l1l2l3); + void CreatePVectorDiagnosticHutter(Vec pg,void* inputs,int analysis_type); + void* GetMatPar(); + + Object* copy(); + void NodeConfiguration(int* beam_node_ids,Node* beam_nodes[2],int* beam_node_offsets); + void MaticeConfiguration(Matice* matice,int matice_offset); + void MatparConfiguration(Matpar* matpar,int matpar_offset); + void ComputePressure(Vec p_g); + void GetNodes(void** vpnodes); + + /*not implemented: */ + int GetShelf(); + int GetOnBed(); + void GetBedList(double*); + void GetThicknessList(double* thickness_list); + void Du(_p_Vec*, double*, double*, void*, int); + void Gradj(_p_Vec*, double*, double*, void*, int, char*); + void GradjDrag(_p_Vec*, double*, double*, void*, int); + void GradjB(_p_Vec*, double*, double*, void*, int); + double Misfit(double*, double*, void*, int); + void GetNodalFunctions(double* l1l2, double gauss_coord); + void GetParameterValue(double* pvalue, double* value_list,double gauss_coord); + void GetJacobianDeterminant(double* pJdet,double* z_list, double gauss_coord); + +}; +#endif /* _BEAM_H */ Index: /issm/trunk/src/c/objects/Matice.cpp =================================================================== --- /issm/trunk/src/c/objects/Matice.cpp (revision 307) +++ /issm/trunk/src/c/objects/Matice.cpp (revision 308) @@ -305,2 +305,11 @@ *pviscosity3d=viscosity3d; } + +double Matice::GetB(){ + return B; +} + +double Matice::GetN(){ + return n; +} + Index: /issm/trunk/src/c/objects/Matice.h =================================================================== --- /issm/trunk/src/c/objects/Matice.h (revision 307) +++ /issm/trunk/src/c/objects/Matice.h (revision 308) @@ -36,4 +36,6 @@ void GetViscosity3d(double* pviscosity3d, double* pepsilon); Object* copy(); + double GetB(); + double GetN(); }; Index: /issm/trunk/src/c/objects/Sing.cpp =================================================================== --- /issm/trunk/src/c/objects/Sing.cpp (revision 308) +++ /issm/trunk/src/c/objects/Sing.cpp (revision 308) @@ -0,0 +1,463 @@ +/*!\file Sing.c + * \brief: implementation of the Sing object + */ + +#ifdef HAVE_CONFIG_H + #include "config.h" +#else +#error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!" +#endif + +#include "stdio.h" +#include "./Sing.h" +#include +#include "../EnumDefinitions/EnumDefinitions.h" +#include "../shared/shared.h" +#include "../DataSet/DataSet.h" +#include "../include/typedefs.h" + + +Sing::Sing(){ + return; +} + +Sing::~Sing(){ + return; +} + +Sing::Sing(int sing_id, int sing_mid, int sing_mparid, int sing_g, double sing_h, double sing_k){ + + id=sing_id; + mid=sing_mid; + mparid=sing_mparid; + + node_id=sing_g; + node_offset=UNDEF; + node=NULL; + + h=sing_h; + k=sing_k; + + matice=NULL; + matice_offset=UNDEF; + matpar=NULL; + matpar_offset=UNDEF; + + return; +} + +#undef __FUNCT__ +#define __FUNCT__ "Sing::Echo" + +void Sing::Echo(void){ + + printf("Sing:\n"); + printf(" id: %i\n",id); + printf(" mid: %i\n",mid); + printf(" mparid: %i\n",mparid); + printf(" node_id=[%i]\n",node_id); + printf(" node_offset=[%i]\n",node_offset); + printf(" matice_offset=%i\n",matice_offset); + printf(" matpar_offset=%i\n",matpar_offset); + printf(" h=[%g]\n",h); + printf(" k=[%g]\n",h); + printf(" node: \n"); + if(node)node->Echo(); + if(matice)matice->Echo(); + if(matpar)matpar->Echo(); + + return; +} + +void Sing::Marshall(char** pmarshalled_dataset){ + + char* marshalled_dataset=NULL; + int enum_type=0; + + /*recover marshalled_dataset: */ + marshalled_dataset=*pmarshalled_dataset; + + /*get enum type of Sing: */ + enum_type=SingEnum(); + + /*marshall enum: */ + memcpy(marshalled_dataset,&enum_type,sizeof(enum_type));marshalled_dataset+=sizeof(enum_type); + + /*marshall Sing data: */ + memcpy(marshalled_dataset,&id,sizeof(id));marshalled_dataset+=sizeof(id); + memcpy(marshalled_dataset,&mid,sizeof(mid));marshalled_dataset+=sizeof(mid); + memcpy(marshalled_dataset,&mparid,sizeof(mparid));marshalled_dataset+=sizeof(mparid); + memcpy(marshalled_dataset,&node_id,sizeof(node_id));marshalled_dataset+=sizeof(node_id); + memcpy(marshalled_dataset,&node,sizeof(node));marshalled_dataset+=sizeof(node); + memcpy(marshalled_dataset,&node_offset,sizeof(node_offset));marshalled_dataset+=sizeof(node_offset); + memcpy(marshalled_dataset,&matice,sizeof(matice));marshalled_dataset+=sizeof(matice); + memcpy(marshalled_dataset,&matice_offset,sizeof(matice_offset));marshalled_dataset+=sizeof(matice_offset); + memcpy(marshalled_dataset,&matpar,sizeof(matpar));marshalled_dataset+=sizeof(matpar); + memcpy(marshalled_dataset,&matpar_offset,sizeof(matpar_offset));marshalled_dataset+=sizeof(matpar_offset); + memcpy(marshalled_dataset,&h,sizeof(h));marshalled_dataset+=sizeof(h); + memcpy(marshalled_dataset,&k,sizeof(k));marshalled_dataset+=sizeof(k); + + *pmarshalled_dataset=marshalled_dataset; + return; +} + +int Sing::MarshallSize(){ + return sizeof(id) + +sizeof(mid) + +sizeof(mparid) + +sizeof(node_id) + +sizeof(node) + +sizeof(node_offset) + +sizeof(matice) + +sizeof(matice_offset) + +sizeof(matpar) + +sizeof(matpar_offset) + +sizeof(h) + +sizeof(k); + +sizeof(int); //sizeof(int) for enum type +} + +char* Sing::GetName(void){ + return "sing"; +} + +void Sing::Demarshall(char** pmarshalled_dataset){ + + char* marshalled_dataset=NULL; + int i; + + /*recover marshalled_dataset: */ + marshalled_dataset=*pmarshalled_dataset; + + /*this time, no need to get enum type, the pointer directly points to the beginning of the + *object data (thanks to DataSet::Demarshall):*/ + + memcpy(&id,marshalled_dataset,sizeof(id));marshalled_dataset+=sizeof(id); + memcpy(&mid,marshalled_dataset,sizeof(mid));marshalled_dataset+=sizeof(mid); + memcpy(&mparid,marshalled_dataset,sizeof(mparid));marshalled_dataset+=sizeof(mparid); + memcpy(&node_id,marshalled_dataset,sizeof(node_id));marshalled_dataset+=sizeof(node_id); + memcpy(&node,marshalled_dataset,sizeof(node));marshalled_dataset+=sizeof(node); + memcpy(&node_offset,marshalled_dataset,sizeof(node_offset));marshalled_dataset+=sizeof(node_offset); + memcpy(&matice,marshalled_dataset,sizeof(matice));marshalled_dataset+=sizeof(matice); + memcpy(&matice_offset,marshalled_dataset,sizeof(matice_offset));marshalled_dataset+=sizeof(matice_offset); + memcpy(&matpar,marshalled_dataset,sizeof(matpar));marshalled_dataset+=sizeof(matpar); + memcpy(&matpar_offset,marshalled_dataset,sizeof(matpar_offset));marshalled_dataset+=sizeof(matpar_offset); + memcpy(&h,marshalled_dataset,sizeof(h));marshalled_dataset+=sizeof(h); + memcpy(&k,marshalled_dataset,sizeof(k));marshalled_dataset+=sizeof(k); + + /*nodes and materials are not pointing to correct objects anymore:*/ + node=NULL; + matice=NULL; + matpar=NULL; + + /*return: */ + *pmarshalled_dataset=marshalled_dataset; + return; +} +int Sing::Enum(void){ + + return SingEnum(); + +} +int Sing::GetId(void){ return id; } + +int Sing::MyRank(void){ + extern int my_rank; + return my_rank; +} + + +#undef __FUNCT__ +#define __FUNCT__ "Sing::Configure" +void Sing::Configure(void* ploadsin, void* pnodesin,void* pmaterialsin){ + + int i; + + DataSet* nodesin=NULL; + DataSet* materialsin=NULL; + + /*Recover pointers :*/ + nodesin=(DataSet*)pnodesin; + materialsin=(DataSet*)pmaterialsin; + + /*Link this element with its nodes, ie find pointers to the nodes in the nodes dataset.: */ + ResolvePointers((Object**)&node,&node_id,&node_offset,1,nodesin); + + /*Same for materials: */ + ResolvePointers((Object**)&matice,&mid,&matice_offset,1,materialsin); + ResolvePointers((Object**)&matpar,&mparid,&matpar_offset,1,materialsin); + +} + +#undef __FUNCT__ +#define __FUNCT__ "Sing::CreateKMatrix" + +void Sing::CreateKMatrix(Mat Kgg,void* inputs,int analysis_type){ + + /*Just branch to the correct element stiffness matrix generator, according to the type of analysis we are carrying out: */ + if (analysis_type==DiagnosticHutterAnalysisEnum()){ + + CreateKMatrixDiagnosticHutter( Kgg,inputs,analysis_type); + + } + else{ + throw ErrorException(__FUNCT__,exprintf("%s%i%s\n","analysis: ",analysis_type," not supported yet")); + } + +} + + +#undef __FUNCT__ +#define __FUNCT__ "Sing::CreateKMatrixDiagnosticHutter" + +void Sing::CreateKMatrixDiagnosticHutter(Mat Kgg,void* vinputs,int analysis_type){ + + const int numgrids=1; + const int NDOF2=2; + const int numdofs=numgrids*NDOF2; + double Ke_gg[numdofs][numdofs]={{1,0},{0,1}}; + int doflist[numdofs]; + int numberofdofspernode; + + ParameterInputs* inputs=NULL; + + /*recover pointers: */ + inputs=(ParameterInputs*)vinputs; + + GetDofList(&doflist[0],&numberofdofspernode); + + MatSetValues(Kgg,numdofs,doflist,numdofs,doflist,(const double*)Ke_gg,ADD_VALUES); + +} + + +#undef __FUNCT__ +#define __FUNCT__ "Sing::CreatePVector" +void Sing::CreatePVector(Vec pg,void* inputs,int analysis_type){ + + /*Just branch to the correct load generator, according to the type of analysis we are carrying out: */ + if (analysis_type==DiagnosticHutterAnalysisEnum()){ + CreatePVectorDiagnosticHutter( pg,inputs,analysis_type); + } + else{ + throw ErrorException(__FUNCT__,exprintf("%s%i%s"," analysis ",analysis_type," not supported yet")); + } + +} + + + +#undef __FUNCT__ +#define __FUNCT__ "Sing::CreatePVectorDiagnosticHutter" + +void Sing::CreatePVectorDiagnosticHutter( Vec pg, void* vinputs, int analysis_type){ + + + const int numgrids=1; + const int NDOF2=2; + const int numdofs=NDOF2*numgrids; + int doflist[numdofs]; + int dofs[1]={0}; + int found=0; + double slope[2]; + double slope2; + double pe_g[numdofs]={0,0}; + double ub,vb; + double constant_part; + int numberofdofspernode; + double rho_ice,gravity,n,B; + + ParameterInputs* inputs=NULL; + + /*recover pointers: */ + inputs=(ParameterInputs*)vinputs; + + found=inputs->Recover("surfaceslopex",&slope[0],1,dofs,numgrids,(void**)&node); + if(!found)throw ErrorException(__FUNCT__," surfaceslopex missing from inputs!"); + found=inputs->Recover("surfaceslopey",&slope[1],1,dofs,numgrids,(void**)&node); + if(!found)throw ErrorException(__FUNCT__," surfaceslopey missing from inputs!"); + + GetDofList(&doflist[0],&numberofdofspernode); + + //compute slope2 + slope2=pow(slope[0],2)+pow(slope[1],2); + + //compute ub + rho_ice=matpar->GetRhoIce(); + gravity=matpar->GetG(); + n=matice->GetN(); + B=matice->GetB(); + + ub=-1.58*pow(10.0,-10.0)*rho_ice*gravity*h*slope[0]; + vb=-1.58*pow(10.0,-10.0)*rho_ice*gravity*h*slope[1]; + + //compute constant_part + constant_part=-2*pow(rho_ice*gravity,n)*pow(slope2,((n-1)/2)); + + pe_g[0]=ub-2.0*pow(rho_ice*gravity,n)*pow(slope2,((n-1)/2.0))*pow(h,n)/(pow(B,n)*(n+1))*slope[0]; + pe_g[1]=vb-2.0*pow(rho_ice*gravity,n)*pow(slope2,((n-1)/2.0))*pow(h,n)/(pow(B,n)*(n+1))*slope[1]; + + VecSetValues(pg,numdofs,doflist,(const double*)pe_g,ADD_VALUES); + + +} + + +#undef __FUNCT__ +#define __FUNCT__ "Sing::UpdateFromInputs" +void Sing::UpdateFromInputs(void* vinputs){ + + int dofs[1]={0}; + double temperature; + double B; + + ParameterInputs* inputs=NULL; + + /*recover pointers: */ + inputs=(ParameterInputs*)vinputs; + + /*Update internal data if inputs holds new values: */ + inputs->Recover("thickness",&h,1,dofs,1,(void**)&node); + inputs->Recover("drag",&k,1,dofs,1,(void**)&node); + + //Update material if necessary + if(inputs->Recover("temperature_average",&temperature,1,dofs,1,(void**)&node)){ + B=Paterson(temperature); + matice->SetB(B); + } + if(inputs->Recover("temperature",&temperature,1,dofs,1,(void**)&node)){ + B=Paterson(temperature); + matice->SetB(B); + } + + if(inputs->Recover("B",&B,1,dofs,1,(void**)&node)){ + matice->SetB(B); + } + +} + +void Sing::GetDofList(int* doflist,int* pnumberofdofspernode){ + + int i; + int doflist_per_node[MAXDOFSPERNODE]; + int numberofdofspernode; + + node->GetDofList(&doflist_per_node[0],&numberofdofspernode); + for(i=0;iGetDofList1(); + +} + + +void* Sing::GetMatPar(){ + return matpar; +} + + +void Sing::GetNodes(void** vpnodes){ + + Node** pnodes=(Node**)vpnodes; + + pnodes[0]=node; +} + + +Object* Sing::copy() { + + return new Sing(*this); + +} + + +#undef __FUNCT__ +#define __FUNCT__ "Sing::NodeConfiguration" +void Sing::NodeConfiguration(int sing_node_id,Node* sing_node,int sing_node_offset){ + + node_id=sing_node_id; + node=sing_node; + node_offset=sing_node_offset; + +} +#undef __FUNCT__ +#define __FUNCT__ "Sing::MaticeConfiguration" +void Sing::MaticeConfiguration(Matice* sing_matice,int sing_matice_offset){ + matice=sing_matice; + matice_offset=sing_matice_offset; +} + +#undef __FUNCT__ +#define __FUNCT__ "Sing::MatparConfiguration" +void Sing::MatparConfiguration(Matpar* sing_matpar,int sing_matpar_offset){ + + matpar=sing_matpar; + matpar_offset=sing_matpar_offset; + +} + + +#undef __FUNCT__ +#define __FUNCT__ "Sing::GetShelf" +int Sing::GetShelf(){ + throw ErrorException(__FUNCT__," not supported yet!"); +} +#undef __FUNCT__ +#define __FUNCT__ "Sing::GetOnBed" +int Sing::GetOnBed(){ + throw ErrorException(__FUNCT__," not supported yet!"); +} + + +#undef __FUNCT__ +#define __FUNCT__ "Sing::GetBedList" +void Sing::GetBedList(double*){ + throw ErrorException(__FUNCT__," not supported yet!"); +} + +#undef __FUNCT__ +#define __FUNCT__ "Sing::Du" +void Sing::Du(_p_Vec*, double*, double*, void*, int){ + throw ErrorException(__FUNCT__," not supported yet!"); +} +#undef __FUNCT__ +#define __FUNCT__ "Sing::Gradj" +void Sing::Gradj(_p_Vec*, double*, double*, void*, int, char*){ + throw ErrorException(__FUNCT__," not supported yet!"); +} +#undef __FUNCT__ +#define __FUNCT__ "Sing::GradjDrag" +void Sing::GradjDrag(_p_Vec*, double*, double*, void*, int){ + throw ErrorException(__FUNCT__," not supported yet!"); +} +#undef __FUNCT__ +#define __FUNCT__ "Sing::GradjB" +void Sing::GradjB(_p_Vec*, double*, double*, void*, int){ + throw ErrorException(__FUNCT__," not supported yet!"); +} +#undef __FUNCT__ +#define __FUNCT__ "Sing::Misfit" +double Sing::Misfit(double*, double*, void*, int){ + throw ErrorException(__FUNCT__," not supported yet!"); +} + +#undef __FUNCT__ +#define __FUNCT__ "Sing::GetThicknessList" +void Sing::GetThicknessList(double* thickness_list){ + throw ErrorException(__FUNCT__," not supported yet!"); +} + +#undef __FUNCT__ +#define __FUNCT__ "Sing::ComputePressure" +void Sing::ComputePressure(Vec p_g){ + throw ErrorException(__FUNCT__," not supported yet!"); +} Index: /issm/trunk/src/c/objects/Sing.h =================================================================== --- /issm/trunk/src/c/objects/Sing.h (revision 308) +++ /issm/trunk/src/c/objects/Sing.h (revision 308) @@ -0,0 +1,85 @@ +/*!\file: Sing.h + * \brief prototypes for Sing element + */ + +#ifndef _SING_H_ +#define _SING_H_ + +#include "./Element.h" +#include "./Node.h" +#include "./Matice.h" +#include "./Matpar.h" +#include "./ParameterInputs.h" + +class Sing: public Element{ + + private: + + /*id:*/ + int id; + + /*node: */ + int node_id; //node id + Node* node; //node pointer + int node_offset; //node offset in nodes dataset + + /*material: */ + int mid; + Matice* matice; + int matice_offset; + + int mparid; + Matpar* matpar; + int matpar_offset; + + double h; + double k; + + public: + + Sing(); + Sing(int sing_id, int sing_mid, int sing_mparid, int sing_g, double sing_h, double sing_k); + ~Sing(); + + void Echo(); + void Marshall(char** pmarshalled_dataset); + int MarshallSize(); + char* GetName(); + void Demarshall(char** pmarshalled_dataset); + int Enum(); + int GetId(); + int MyRank(); + void Configure(void* loads,void* nodes,void* materials); + void CreateKMatrix(Mat Kgg,void* inputs,int analysis_type); + void CreatePVector(Vec pg, void* inputs, int analysis_type); + void UpdateFromInputs(void* inputs); + void GetDofList(int* doflist,int* pnumberofdofs); + void GetDofList1(int* doflist); + void CreateKMatrixDiagnosticHutter(Mat Kgg,void* inputs,int analysis_type); + void GetParameterValue(double* pp, double* plist, double* gauss_l1l2l3); + void CreatePVectorDiagnosticHutter(Vec pg,void* inputs,int analysis_type); + void* GetMatPar(); + + Object* copy(); + void NodeConfiguration(int sing_node_id,Node* sing_node,int sing_node_offset); + void MaticeConfiguration(Matice* matice,int matice_offset); + void MatparConfiguration(Matpar* matpar,int matpar_offset); + void ComputePressure(Vec p_g); + void GetNodes(void** vpnodes); + + /*not implemented: */ + int GetShelf(); + int GetOnBed(); + void GetBedList(double*); + void GetThicknessList(double* thickness_list); + void Du(_p_Vec*, double*, double*, void*, int); + void Gradj(_p_Vec*, double*, double*, void*, int, char*); + void GradjDrag(_p_Vec*, double*, double*, void*, int); + void GradjB(_p_Vec*, double*, double*, void*, int); + double Misfit(double*, double*, void*, int); + + + +}; +#endif //ifndef _SING_H_ + Index: /issm/trunk/src/c/objects/objects.h =================================================================== --- /issm/trunk/src/c/objects/objects.h (revision 307) +++ /issm/trunk/src/c/objects/objects.h (revision 308) @@ -16,4 +16,6 @@ #include "./Penta.h" #include "./Tria.h" +#include "./Sing.h" +#include "./Beam.h" #include "./Spc.h" #include "./Rgb.h" Index: /issm/trunk/src/m/classes/public/marshall.m =================================================================== --- /issm/trunk/src/m/classes/public/marshall.m (revision 307) +++ /issm/trunk/src/m/classes/public/marshall.m (revision 308) @@ -32,4 +32,5 @@ WriteData(fid,md.elements,'Mat','elements'); WriteData(fid,md.elements_type,'Mat','elements_type'); +WriteData(fid,md.gridonhutter,'Mat','gridonhutter'); if strcmpi(md.type,'3d'),