Index: sm/trunk/src/c/modules/ModelProcessorx/Balancedthickness/CreateElementsNodesAndMaterialsBalancedthickness.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Balancedthickness/CreateElementsNodesAndMaterialsBalancedthickness.cpp (revision 3997) +++ (revision ) @@ -1,186 +1,0 @@ -/* - * CreateElementsNodesAndMaterialsBalancedthickness.c: - */ - -#include "../../../DataSet/DataSet.h" -#include "../../../toolkits/toolkits.h" -#include "../../../io/io.h" -#include "../../../EnumDefinitions/EnumDefinitions.h" -#include "../../../objects/objects.h" -#include "../../../shared/shared.h" -#include "../../../include/include.h" -#include "../../MeshPartitionx/MeshPartitionx.h" -#include "../ModelProcessorx.h" - -void CreateElementsNodesAndMaterialsBalancedthickness(DataSet** pelements,DataSet** pnodes, DataSet** pvertices, DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ - - /*Intermediary*/ - int i,j,k,n; - - /*DataSets: */ - DataSet* elements = NULL; - DataSet* nodes = NULL; - DataSet* vertices = NULL; - DataSet* materials = NULL; - - /*First create the elements, nodes and material properties: */ - elements = new DataSet(ElementsEnum); - nodes = new DataSet(NodesEnum); - vertices = new DataSet(VerticesEnum); - materials = new DataSet(MaterialsEnum); - - /*Partition elements and vertices and nodes: */ - Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); - - /*2d mesh: */ - if (iomodel->dim==2){ - - /*Fetch data needed: */ - IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); - IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); - IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); - IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); - IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate"); - IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); - IoModelFetchData(&iomodel->dhdt,NULL,NULL,iomodel_handle,"dhdt"); - - for (i=0;inumberofelements;i++){ - - if(iomodel->my_elements[i]){ - - /*Create and add tria element to elements dataset: */ - elements->AddObject(new Tria(i+1,i,iomodel)); - - /*Create and add material property to materials dataset: */ - materials->AddObject(new Matice(i+1,i,iomodel,3)); - } - }//for (i=0;ielements); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->surface); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->elementoniceshelf); - xfree((void**)&iomodel->elementonwater); - xfree((void**)&iomodel->vx); - xfree((void**)&iomodel->vy); - xfree((void**)&iomodel->melting_rate); - xfree((void**)&iomodel->accumulation_rate); - xfree((void**)&iomodel->dhdt); - - } - else{ // if (dim==2) - - /*Fetch data needed: */ - IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); - IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed"); - IoModelFetchData(&iomodel->elementonsurface,NULL,NULL,iomodel_handle,"elementonsurface"); - IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); - IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); - IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); - IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate"); - IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); - IoModelFetchData(&iomodel->dhdt,NULL,NULL,iomodel_handle,"dhdt"); - IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements"); - IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements"); - - for (i=0;inumberofelements;i++){ - if(iomodel->my_elements[i]){ - /*Create and add penta element to elements dataset: */ - elements->AddObject(new Penta(i+1,i,iomodel)); - - /*Create and add material property to materials dataset: */ - materials->AddObject(new Matice(i+1,i,iomodel,6)); - } - }//for (i=0;ielements); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->surface); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->elementoniceshelf); - xfree((void**)&iomodel->elementonbed); - xfree((void**)&iomodel->elementonsurface); - xfree((void**)&iomodel->elementonwater); - xfree((void**)&iomodel->vx); - xfree((void**)&iomodel->vy); - xfree((void**)&iomodel->melting_rate); - xfree((void**)&iomodel->accumulation_rate); - xfree((void**)&iomodel->dhdt); - xfree((void**)&iomodel->upperelements); - xfree((void**)&iomodel->lowerelements); - - - } //if (dim==2) - - /*Add new constrant material property to materials, at the end: */ - materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel)); - - /*First fetch data: */ - if (iomodel->dim==3){ - IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); - IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); - } - IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); - IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); - IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); - IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); - IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); - IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); - - for (i=0;inumberofvertices;i++){ - - /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ - if(iomodel->my_vertices[i]){ - - /*Add vertex to vertices dataset: */ - vertices->AddObject(new Vertex(i+1,i,iomodel)); - - /*Add node to nodes dataset: */ - nodes->AddObject(new Node(i+1,i,iomodel)); - - } - } - - /*Clean fetched data: */ - xfree((void**)&iomodel->deadgrids); - xfree((void**)&iomodel->x); - xfree((void**)&iomodel->y); - xfree((void**)&iomodel->z); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->gridonbed); - xfree((void**)&iomodel->gridonsurface); - xfree((void**)&iomodel->uppernodes); - xfree((void**)&iomodel->gridonicesheet); - xfree((void**)&iomodel->gridoniceshelf); - - /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these - * datasets, it will not be redone: */ - elements->Presort(); - vertices->Presort(); - nodes->Presort(); - materials->Presort(); - - cleanup_and_return: - - /*Assign output pointer: */ - *pelements=elements; - *pnodes=nodes; - *pvertices=vertices; - *pmaterials=materials; - -} Index: /issm/trunk/src/c/modules/ModelProcessorx/Balancedthickness/CreateNodesBalancedthickness.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Balancedthickness/CreateNodesBalancedthickness.cpp (revision 3998) +++ /issm/trunk/src/c/modules/ModelProcessorx/Balancedthickness/CreateNodesBalancedthickness.cpp (revision 3998) @@ -0,0 +1,186 @@ +/* + * CreateElementsNodesAndMaterialsBalancedthickness.c: + */ + +#include "../../../DataSet/DataSet.h" +#include "../../../toolkits/toolkits.h" +#include "../../../io/io.h" +#include "../../../EnumDefinitions/EnumDefinitions.h" +#include "../../../objects/objects.h" +#include "../../../shared/shared.h" +#include "../../../include/include.h" +#include "../../MeshPartitionx/MeshPartitionx.h" +#include "../ModelProcessorx.h" + +void CreateElementsNodesAndMaterialsBalancedthickness(DataSet** pelements,DataSet** pnodes, DataSet** pvertices, DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ + + /*Intermediary*/ + int i,j,k,n; + + /*DataSets: */ + DataSet* elements = NULL; + DataSet* nodes = NULL; + DataSet* vertices = NULL; + DataSet* materials = NULL; + + /*First create the elements, nodes and material properties: */ + elements = new DataSet(ElementsEnum); + nodes = new DataSet(NodesEnum); + vertices = new DataSet(VerticesEnum); + materials = new DataSet(MaterialsEnum); + + /*Partition elements and vertices and nodes: */ + Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); + + /*2d mesh: */ + if (iomodel->dim==2){ + + /*Fetch data needed: */ + IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); + IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); + IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); + IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); + IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate"); + IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); + IoModelFetchData(&iomodel->dhdt,NULL,NULL,iomodel_handle,"dhdt"); + + for (i=0;inumberofelements;i++){ + + if(iomodel->my_elements[i]){ + + /*Create and add tria element to elements dataset: */ + elements->AddObject(new Tria(i+1,i,iomodel)); + + /*Create and add material property to materials dataset: */ + materials->AddObject(new Matice(i+1,i,iomodel,3)); + } + }//for (i=0;ielements); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->surface); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->elementoniceshelf); + xfree((void**)&iomodel->elementonwater); + xfree((void**)&iomodel->vx); + xfree((void**)&iomodel->vy); + xfree((void**)&iomodel->melting_rate); + xfree((void**)&iomodel->accumulation_rate); + xfree((void**)&iomodel->dhdt); + + } + else{ // if (dim==2) + + /*Fetch data needed: */ + IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); + IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed"); + IoModelFetchData(&iomodel->elementonsurface,NULL,NULL,iomodel_handle,"elementonsurface"); + IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); + IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); + IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); + IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate"); + IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); + IoModelFetchData(&iomodel->dhdt,NULL,NULL,iomodel_handle,"dhdt"); + IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements"); + IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements"); + + for (i=0;inumberofelements;i++){ + if(iomodel->my_elements[i]){ + /*Create and add penta element to elements dataset: */ + elements->AddObject(new Penta(i+1,i,iomodel)); + + /*Create and add material property to materials dataset: */ + materials->AddObject(new Matice(i+1,i,iomodel,6)); + } + }//for (i=0;ielements); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->surface); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->elementoniceshelf); + xfree((void**)&iomodel->elementonbed); + xfree((void**)&iomodel->elementonsurface); + xfree((void**)&iomodel->elementonwater); + xfree((void**)&iomodel->vx); + xfree((void**)&iomodel->vy); + xfree((void**)&iomodel->melting_rate); + xfree((void**)&iomodel->accumulation_rate); + xfree((void**)&iomodel->dhdt); + xfree((void**)&iomodel->upperelements); + xfree((void**)&iomodel->lowerelements); + + + } //if (dim==2) + + /*Add new constrant material property to materials, at the end: */ + materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel)); + + /*First fetch data: */ + if (iomodel->dim==3){ + IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); + IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); + } + IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); + IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); + IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); + IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); + IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); + IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); + + for (i=0;inumberofvertices;i++){ + + /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ + if(iomodel->my_vertices[i]){ + + /*Add vertex to vertices dataset: */ + vertices->AddObject(new Vertex(i+1,i,iomodel)); + + /*Add node to nodes dataset: */ + nodes->AddObject(new Node(i+1,i,iomodel)); + + } + } + + /*Clean fetched data: */ + xfree((void**)&iomodel->deadgrids); + xfree((void**)&iomodel->x); + xfree((void**)&iomodel->y); + xfree((void**)&iomodel->z); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->gridonbed); + xfree((void**)&iomodel->gridonsurface); + xfree((void**)&iomodel->uppernodes); + xfree((void**)&iomodel->gridonicesheet); + xfree((void**)&iomodel->gridoniceshelf); + + /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these + * datasets, it will not be redone: */ + elements->Presort(); + vertices->Presort(); + nodes->Presort(); + materials->Presort(); + + cleanup_and_return: + + /*Assign output pointer: */ + *pelements=elements; + *pnodes=nodes; + *pvertices=vertices; + *pmaterials=materials; + +} Index: sm/trunk/src/c/modules/ModelProcessorx/Balancedthickness2/CreateElementsNodesAndMaterialsBalancedthickness2.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Balancedthickness2/CreateElementsNodesAndMaterialsBalancedthickness2.cpp (revision 3997) +++ (revision ) @@ -1,163 +1,0 @@ -/* - * CreateElementsNodesAndMaterialsBalancedthickness2.c: - */ - -#include "../../../DataSet/DataSet.h" -#include "../../../toolkits/toolkits.h" -#include "../../../io/io.h" -#include "../../../EnumDefinitions/EnumDefinitions.h" -#include "../../../objects/objects.h" -#include "../../../shared/shared.h" -#include "../../../include/include.h" -#include "../../MeshPartitionx/MeshPartitionx.h" -#include "../ModelProcessorx.h" - -void CreateElementsNodesAndMaterialsBalancedthickness2(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ - - /*Intermediary*/ - int i,j; - int vertex_index; - int node_index; - - /*DataSets: */ - DataSet* elements = NULL; - DataSet* nodes = NULL; - DataSet* vertices = NULL; - DataSet* materials = NULL; - - /*First create the elements, nodes and material properties: */ - elements = new DataSet(ElementsEnum); - nodes = new DataSet(NodesEnum); - vertices = new DataSet(VerticesEnum); - materials = new DataSet(MaterialsEnum); - - /*Partition elements and vertices and nodes: */ - Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); - - /*elements created vary if we are dealing with a 2d mesh, or a 3d mesh: */ - /*2d mesh: */ - if (iomodel->dim==2){ - - /*Fetch data needed: */ - IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); - IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); - IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); - IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); - IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"melting_rate"); - IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); - IoModelFetchData(&iomodel->dhdt,NULL,NULL,iomodel_handle,"dhdt"); - - for (i=0;inumberofelements;i++){ - - if(iomodel->my_elements[i]){ - - /*Create and add tria element to elements dataset: */ - elements->AddObject(new Tria(i+1,i,iomodel)); - - /*Create and add material property to materials dataset: */ - materials->AddObject(new Matice(i+1,i,iomodel,3)); - } - }//for (i=0;ithickness); - xfree((void**)&iomodel->surface); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->elementoniceshelf); - xfree((void**)&iomodel->elementonwater); - xfree((void**)&iomodel->vx); - xfree((void**)&iomodel->vy); - xfree((void**)&iomodel->melting_rate); - xfree((void**)&iomodel->accumulation_rate); - xfree((void**)&iomodel->dhdt); - - } - else{ // if (dim==2) - ISSMERROR("not implemented yet"); - } //if (dim==2) - - /*Add new constrant material property tgo materials, at the end: */ - materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel)); - - /*Create nodes and vertices: */ - IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); - IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); - IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); - IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); - IoModelFetchData(&iomodel->gridonhutter,NULL,NULL,iomodel_handle,"gridonhutter"); - IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); - IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); - if (iomodel->dim==3){ - IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); - IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); - } - - /*Build Vertices dataset*/ - for (i=0;inumberofvertices;i++){ - - /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ - if(iomodel->my_vertices[i]){ - - /*Add vertex to vertices dataset: */ - vertices->AddObject(new Vertex(i+1,i,iomodel)); - - } - } - - /*Build Nodes dataset -> 3 for each element (Discontinuous Galerkin)*/ - for (i=0;inumberofelements;i++){ - for (j=0;j<3;j++){ - - if(iomodel->my_nodes[3*i+j]){ - - //Get index of the vertex on which the current node is located - vertex_index=(int)*(iomodel->elements+3*i+j)-1; //(Matlab to C indexing) - ISSMASSERT(vertex_index>=0 && vertex_indexnumberofvertices); - - //Compute Node index (id-1) - node_index=3*i+j; - - /*Add node to nodes dataset: */ - nodes->AddObject(new Node(i+1,vertex_index,node_index,iomodel)); - - } - } - } - - /*Clean fetched data: */ - xfree((void**)&iomodel->deadgrids); - xfree((void**)&iomodel->x); - xfree((void**)&iomodel->y); - xfree((void**)&iomodel->z); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->gridonbed); - xfree((void**)&iomodel->gridonsurface); - xfree((void**)&iomodel->gridonhutter); - xfree((void**)&iomodel->uppernodes); - xfree((void**)&iomodel->gridonicesheet); - xfree((void**)&iomodel->gridoniceshelf); - xfree((void**)&iomodel->elements); - - /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these - * datasets, it will not be redone: */ - elements->Presort(); - nodes->Presort(); - vertices->Presort(); - materials->Presort(); - - cleanup_and_return: - - /*Assign output pointer: */ - *pelements=elements; - *pnodes=nodes; - *pvertices=vertices; - *pmaterials=materials; -} Index: /issm/trunk/src/c/modules/ModelProcessorx/Balancedthickness2/CreateNodesBalancedthickness2.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Balancedthickness2/CreateNodesBalancedthickness2.cpp (revision 3998) +++ /issm/trunk/src/c/modules/ModelProcessorx/Balancedthickness2/CreateNodesBalancedthickness2.cpp (revision 3998) @@ -0,0 +1,163 @@ +/* + * CreateElementsNodesAndMaterialsBalancedthickness2.c: + */ + +#include "../../../DataSet/DataSet.h" +#include "../../../toolkits/toolkits.h" +#include "../../../io/io.h" +#include "../../../EnumDefinitions/EnumDefinitions.h" +#include "../../../objects/objects.h" +#include "../../../shared/shared.h" +#include "../../../include/include.h" +#include "../../MeshPartitionx/MeshPartitionx.h" +#include "../ModelProcessorx.h" + +void CreateElementsNodesAndMaterialsBalancedthickness2(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ + + /*Intermediary*/ + int i,j; + int vertex_index; + int node_index; + + /*DataSets: */ + DataSet* elements = NULL; + DataSet* nodes = NULL; + DataSet* vertices = NULL; + DataSet* materials = NULL; + + /*First create the elements, nodes and material properties: */ + elements = new DataSet(ElementsEnum); + nodes = new DataSet(NodesEnum); + vertices = new DataSet(VerticesEnum); + materials = new DataSet(MaterialsEnum); + + /*Partition elements and vertices and nodes: */ + Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); + + /*elements created vary if we are dealing with a 2d mesh, or a 3d mesh: */ + /*2d mesh: */ + if (iomodel->dim==2){ + + /*Fetch data needed: */ + IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); + IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); + IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); + IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); + IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"melting_rate"); + IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); + IoModelFetchData(&iomodel->dhdt,NULL,NULL,iomodel_handle,"dhdt"); + + for (i=0;inumberofelements;i++){ + + if(iomodel->my_elements[i]){ + + /*Create and add tria element to elements dataset: */ + elements->AddObject(new Tria(i+1,i,iomodel)); + + /*Create and add material property to materials dataset: */ + materials->AddObject(new Matice(i+1,i,iomodel,3)); + } + }//for (i=0;ithickness); + xfree((void**)&iomodel->surface); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->elementoniceshelf); + xfree((void**)&iomodel->elementonwater); + xfree((void**)&iomodel->vx); + xfree((void**)&iomodel->vy); + xfree((void**)&iomodel->melting_rate); + xfree((void**)&iomodel->accumulation_rate); + xfree((void**)&iomodel->dhdt); + + } + else{ // if (dim==2) + ISSMERROR("not implemented yet"); + } //if (dim==2) + + /*Add new constrant material property tgo materials, at the end: */ + materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel)); + + /*Create nodes and vertices: */ + IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); + IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); + IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); + IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); + IoModelFetchData(&iomodel->gridonhutter,NULL,NULL,iomodel_handle,"gridonhutter"); + IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); + IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); + if (iomodel->dim==3){ + IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); + IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); + } + + /*Build Vertices dataset*/ + for (i=0;inumberofvertices;i++){ + + /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ + if(iomodel->my_vertices[i]){ + + /*Add vertex to vertices dataset: */ + vertices->AddObject(new Vertex(i+1,i,iomodel)); + + } + } + + /*Build Nodes dataset -> 3 for each element (Discontinuous Galerkin)*/ + for (i=0;inumberofelements;i++){ + for (j=0;j<3;j++){ + + if(iomodel->my_nodes[3*i+j]){ + + //Get index of the vertex on which the current node is located + vertex_index=(int)*(iomodel->elements+3*i+j)-1; //(Matlab to C indexing) + ISSMASSERT(vertex_index>=0 && vertex_indexnumberofvertices); + + //Compute Node index (id-1) + node_index=3*i+j; + + /*Add node to nodes dataset: */ + nodes->AddObject(new Node(i+1,vertex_index,node_index,iomodel)); + + } + } + } + + /*Clean fetched data: */ + xfree((void**)&iomodel->deadgrids); + xfree((void**)&iomodel->x); + xfree((void**)&iomodel->y); + xfree((void**)&iomodel->z); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->gridonbed); + xfree((void**)&iomodel->gridonsurface); + xfree((void**)&iomodel->gridonhutter); + xfree((void**)&iomodel->uppernodes); + xfree((void**)&iomodel->gridonicesheet); + xfree((void**)&iomodel->gridoniceshelf); + xfree((void**)&iomodel->elements); + + /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these + * datasets, it will not be redone: */ + elements->Presort(); + nodes->Presort(); + vertices->Presort(); + materials->Presort(); + + cleanup_and_return: + + /*Assign output pointer: */ + *pelements=elements; + *pnodes=nodes; + *pvertices=vertices; + *pmaterials=materials; +} Index: sm/trunk/src/c/modules/ModelProcessorx/Balancedvelocities/CreateElementsNodesAndMaterialsBalancedvelocities.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Balancedvelocities/CreateElementsNodesAndMaterialsBalancedvelocities.cpp (revision 3997) +++ (revision ) @@ -1,186 +1,0 @@ -/* - * CreateElementsNodesAndMaterialsBalancedvelocities.c: - */ - -#include "../../../DataSet/DataSet.h" -#include "../../../toolkits/toolkits.h" -#include "../../../io/io.h" -#include "../../../EnumDefinitions/EnumDefinitions.h" -#include "../../../objects/objects.h" -#include "../../../shared/shared.h" -#include "../../MeshPartitionx/MeshPartitionx.h" -#include "../../../include/include.h" -#include "../ModelProcessorx.h" - -void CreateElementsNodesAndMaterialsBalancedvelocities(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ - - /*Intermediary*/ - int i,j,k,n; - - /*DataSets: */ - DataSet* elements = NULL; - DataSet* nodes = NULL; - DataSet* vertices = NULL; - DataSet* materials = NULL; - - /*First create the elements, nodes and material properties: */ - elements = new DataSet(ElementsEnum); - nodes = new DataSet(NodesEnum); - vertices = new DataSet(VerticesEnum); - materials = new DataSet(MaterialsEnum); - - /*Partition elements and vertices and nodes: */ - Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); - - /*2d mesh: */ - if (iomodel->dim==2){ - - /*Fetch data needed: */ - IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); - IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); - IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); - IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); - IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"melting_rate"); - IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"accumulation_rate"); - - for (i=0;inumberofelements;i++){ - - if(iomodel->my_elements[i]){ - - /*Create and add tria element to elements dataset: */ - elements->AddObject(new Tria(i+1,i,iomodel)); - - /*Create and add material property to materials dataset: */ - materials->AddObject(new Matice(i+1,i,iomodel,3)); - } - - }//for (i=0;ielements); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->surface); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->elementoniceshelf); - xfree((void**)&iomodel->elementonwater); - xfree((void**)&iomodel->vx); - xfree((void**)&iomodel->vy); - xfree((void**)&iomodel->melting_rate); - xfree((void**)&iomodel->accumulation_rate); - - } - else{ // if (dim==2) - - /*Fetch data needed: */ - IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); - IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed"); - IoModelFetchData(&iomodel->elementonsurface,NULL,NULL,iomodel_handle,"elementonsurface"); - IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); - IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); - IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); - IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate"); - IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); - IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements"); - IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements"); - - for (i=0;inumberofelements;i++){ - if(iomodel->my_elements[i]){ - /*Create and add penta element to elements dataset: */ - elements->AddObject(new Penta(i+1,i,iomodel)); - - /*Create and add material property to materials dataset: */ - materials->AddObject(new Matice(i+1,i,iomodel,6)); - - }//if(my_elements[i]) - - }//for (i=0;ielements); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->surface); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->elementoniceshelf); - xfree((void**)&iomodel->elementonbed); - xfree((void**)&iomodel->elementonsurface); - xfree((void**)&iomodel->elementonwater); - xfree((void**)&iomodel->vx); - xfree((void**)&iomodel->vy); - xfree((void**)&iomodel->melting_rate); - xfree((void**)&iomodel->accumulation_rate); - xfree((void**)&iomodel->upperelements); - xfree((void**)&iomodel->lowerelements); - - - } //if (dim==2) - - /*Add new constrant material property to materials, at the end: */ - materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel)); //put it at the end of the materials - - /*First fetch data: */ - if (iomodel->dim==3){ - IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); - IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); - } - IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); - IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); - IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); - IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); - IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); - IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); - - for (i=0;inumberofvertices;i++){ - - /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ - if(iomodel->my_vertices[i]){ - - /*Add vertex to vertices dataset: */ - vertices->AddObject(new Vertex(i+1,i,iomodel)); - - /*Add node to nodes dataset: */ - nodes->AddObject(new Node(i+1,i,iomodel)); - - } - } - - /*Clean fetched data: */ - xfree((void**)&iomodel->deadgrids); - xfree((void**)&iomodel->x); - xfree((void**)&iomodel->y); - xfree((void**)&iomodel->z); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->gridonbed); - xfree((void**)&iomodel->gridonsurface); - xfree((void**)&iomodel->uppernodes); - xfree((void**)&iomodel->gridonicesheet); - xfree((void**)&iomodel->gridoniceshelf); - - /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these - * datasets, it will not be redone: */ - elements->Presort(); - nodes->Presort(); - vertices->Presort(); - materials->Presort(); - - cleanup_and_return: - - /*Assign output pointer: */ - *pelements=elements; - *pnodes=nodes; - *pvertices=vertices; - *pmaterials=materials; - -} Index: /issm/trunk/src/c/modules/ModelProcessorx/Balancedvelocities/CreateNodesBalancedvelocities.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Balancedvelocities/CreateNodesBalancedvelocities.cpp (revision 3998) +++ /issm/trunk/src/c/modules/ModelProcessorx/Balancedvelocities/CreateNodesBalancedvelocities.cpp (revision 3998) @@ -0,0 +1,186 @@ +/* + * CreateElementsNodesAndMaterialsBalancedvelocities.c: + */ + +#include "../../../DataSet/DataSet.h" +#include "../../../toolkits/toolkits.h" +#include "../../../io/io.h" +#include "../../../EnumDefinitions/EnumDefinitions.h" +#include "../../../objects/objects.h" +#include "../../../shared/shared.h" +#include "../../MeshPartitionx/MeshPartitionx.h" +#include "../../../include/include.h" +#include "../ModelProcessorx.h" + +void CreateElementsNodesAndMaterialsBalancedvelocities(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ + + /*Intermediary*/ + int i,j,k,n; + + /*DataSets: */ + DataSet* elements = NULL; + DataSet* nodes = NULL; + DataSet* vertices = NULL; + DataSet* materials = NULL; + + /*First create the elements, nodes and material properties: */ + elements = new DataSet(ElementsEnum); + nodes = new DataSet(NodesEnum); + vertices = new DataSet(VerticesEnum); + materials = new DataSet(MaterialsEnum); + + /*Partition elements and vertices and nodes: */ + Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); + + /*2d mesh: */ + if (iomodel->dim==2){ + + /*Fetch data needed: */ + IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); + IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); + IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); + IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); + IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"melting_rate"); + IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"accumulation_rate"); + + for (i=0;inumberofelements;i++){ + + if(iomodel->my_elements[i]){ + + /*Create and add tria element to elements dataset: */ + elements->AddObject(new Tria(i+1,i,iomodel)); + + /*Create and add material property to materials dataset: */ + materials->AddObject(new Matice(i+1,i,iomodel,3)); + } + + }//for (i=0;ielements); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->surface); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->elementoniceshelf); + xfree((void**)&iomodel->elementonwater); + xfree((void**)&iomodel->vx); + xfree((void**)&iomodel->vy); + xfree((void**)&iomodel->melting_rate); + xfree((void**)&iomodel->accumulation_rate); + + } + else{ // if (dim==2) + + /*Fetch data needed: */ + IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); + IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed"); + IoModelFetchData(&iomodel->elementonsurface,NULL,NULL,iomodel_handle,"elementonsurface"); + IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); + IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); + IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); + IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate"); + IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); + IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements"); + IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements"); + + for (i=0;inumberofelements;i++){ + if(iomodel->my_elements[i]){ + /*Create and add penta element to elements dataset: */ + elements->AddObject(new Penta(i+1,i,iomodel)); + + /*Create and add material property to materials dataset: */ + materials->AddObject(new Matice(i+1,i,iomodel,6)); + + }//if(my_elements[i]) + + }//for (i=0;ielements); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->surface); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->elementoniceshelf); + xfree((void**)&iomodel->elementonbed); + xfree((void**)&iomodel->elementonsurface); + xfree((void**)&iomodel->elementonwater); + xfree((void**)&iomodel->vx); + xfree((void**)&iomodel->vy); + xfree((void**)&iomodel->melting_rate); + xfree((void**)&iomodel->accumulation_rate); + xfree((void**)&iomodel->upperelements); + xfree((void**)&iomodel->lowerelements); + + + } //if (dim==2) + + /*Add new constrant material property to materials, at the end: */ + materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel)); //put it at the end of the materials + + /*First fetch data: */ + if (iomodel->dim==3){ + IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); + IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); + } + IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); + IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); + IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); + IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); + IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); + IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); + + for (i=0;inumberofvertices;i++){ + + /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ + if(iomodel->my_vertices[i]){ + + /*Add vertex to vertices dataset: */ + vertices->AddObject(new Vertex(i+1,i,iomodel)); + + /*Add node to nodes dataset: */ + nodes->AddObject(new Node(i+1,i,iomodel)); + + } + } + + /*Clean fetched data: */ + xfree((void**)&iomodel->deadgrids); + xfree((void**)&iomodel->x); + xfree((void**)&iomodel->y); + xfree((void**)&iomodel->z); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->gridonbed); + xfree((void**)&iomodel->gridonsurface); + xfree((void**)&iomodel->uppernodes); + xfree((void**)&iomodel->gridonicesheet); + xfree((void**)&iomodel->gridoniceshelf); + + /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these + * datasets, it will not be redone: */ + elements->Presort(); + nodes->Presort(); + vertices->Presort(); + materials->Presort(); + + cleanup_and_return: + + /*Assign output pointer: */ + *pelements=elements; + *pnodes=nodes; + *pvertices=vertices; + *pmaterials=materials; + +} Index: /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHoriz/CreateNodesDiagnosticHoriz.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHoriz/CreateNodesDiagnosticHoriz.cpp (revision 3997) +++ /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHoriz/CreateNodesDiagnosticHoriz.cpp (revision 3998) @@ -17,5 +17,5 @@ /*Intermediary*/ int i; - bool galerkin_continuous=true; + bool continuous_galerkin=true; /*DataSets: */ @@ -31,6 +31,6 @@ totalnodes=iomodel->nodecounter; - /*Galerkin continus partition of nodes: */ - NodesPartitioning(&iomodel->my_nodes,iomodel->my_elements, iomodel->my_vertices, iomodel->my_bordervertices, iomodel, iomodel_handle,galerkin_continuous); + /*Continuous Galerkin partition of nodes: */ + NodesPartitioning(&iomodel->my_nodes,iomodel->my_elements, iomodel->my_vertices, iomodel->my_bordervertices, iomodel, iomodel_handle,continuous_galerkin); /*Create nodes: */ Index: sm/trunk/src/c/modules/ModelProcessorx/DiagnosticHutter/CreateElementsNodesAndMaterialsDiagnosticHutter.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHutter/CreateElementsNodesAndMaterialsDiagnosticHutter.cpp (revision 3997) +++ (revision ) @@ -1,188 +1,0 @@ -/* - * CreateElementsNodesAndMaterialsDiagnosticHutter.c: - */ - -#include "../../../DataSet/DataSet.h" -#include "../../../toolkits/toolkits.h" -#include "../../../io/io.h" -#include "../../../EnumDefinitions/EnumDefinitions.h" -#include "../../../objects/objects.h" -#include "../../../shared/shared.h" -#include "../../MeshPartitionx/MeshPartitionx.h" -#include "../../../include/include.h" -#include "../ModelProcessorx.h" - -void CreateElementsNodesAndMaterialsDiagnosticHutter(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ - - /*Intermediary*/ - int i,j,k,n; - - /*DataSets: */ - DataSet* elements = NULL; - DataSet* nodes = NULL; - DataSet* vertices = NULL; - DataSet* materials = NULL; - - /*First create the elements, nodes and material properties: */ - elements = new DataSet(ElementsEnum); - nodes = new DataSet(NodesEnum); - vertices = new DataSet(VerticesEnum); - materials = new DataSet(MaterialsEnum); - - /*Now, is the flag ishutter on? otherwise, do nothing: */ - if (!iomodel->ishutter)goto cleanup_and_return; - - /*Partition elements and vertices and nodes: */ - Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); - - /*Hutter elements can be partitioned using epart, even if - * each hutter elements either lies on a node (in 2d), or a pair of vertically juxtaposed nodes (in 3d): */ - - /*Fetch data temporarily needed: */ - IoModelFetchData(&iomodel->gridonhutter,NULL,NULL,iomodel_handle,"gridonhutter"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); - IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); - IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); - IoModelFetchData(&iomodel->drag_coefficient,NULL,NULL,iomodel_handle,"drag_coefficient"); - IoModelFetchData(&iomodel->rheology_B,NULL,NULL,iomodel_handle,"rheology_B"); - IoModelFetchData(&iomodel->rheology_n,NULL,NULL,iomodel_handle,"rheology_n"); - IoModelFetchData(&iomodel->elements_type,NULL,NULL,iomodel_handle,"elements_type"); - IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed"); - IoModelFetchData(&iomodel->elementonsurface,NULL,NULL,iomodel_handle,"elementonsurface"); - IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); - IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); - IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements"); - IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements"); - - /*2d mesh: */ - if (iomodel->dim==2){ - - for (i=0;inumberofelements;i++){ - - if(iomodel->my_elements[i]){ - - if (*(iomodel->elements_type+2*i+0)==HutterFormulationEnum){ //create only Hutter elements - - /*Create and add penta element to elements dataset: */ - elements->AddObject(new Tria(i+1,i,iomodel)); - - /*Create and add material property to materials dataset: */ - materials->AddObject(new Matice(i+1,i,iomodel,3)); - } - } - } //for (i=0;inumberofvertices;i++) - } //if (iomodel->dim==2) - else{ - - for (i=0;inumberofelements;i++){ - - if(iomodel->my_elements[i]){ - - if (*(iomodel->elements_type+2*i+0)==HutterFormulationEnum){ //create only Hutter elements - - /*Create and add penta element to elements dataset: */ - elements->AddObject(new Penta(i+1,i,iomodel)); - - - /*Create and add material property to materials dataset: */ - materials->AddObject(new Matice(i+1,i,iomodel,6)); - - } - } - - } //for (i=0;inumberofelements;i++) - - } //if (iomodel->dim==2) - - /*Free data: */ - xfree((void**)&iomodel->elements); - xfree((void**)&iomodel->elementonbed); - xfree((void**)&iomodel->elementonsurface); - xfree((void**)&iomodel->elements_type); - xfree((void**)&iomodel->gridonhutter); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->surface); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->gridonbed); - xfree((void**)&iomodel->elementonwater); - xfree((void**)&iomodel->uppernodes); - xfree((void**)&iomodel->drag_coefficient); - xfree((void**)&iomodel->rheology_B); - xfree((void**)&iomodel->rheology_n); - xfree((void**)&iomodel->upperelements); - xfree((void**)&iomodel->lowerelements); - - /*Add new constrant material property to materials, at the end: */ - if (iomodel->dim==2){ - materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel)); //put it at the end of the materials - } - else{ - materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel)); //put it at the end of the materials - } - - /*First fetch data: */ - if (iomodel->dim==3){ - IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); - IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); - } - IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); - IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); - IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); - IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); - IoModelFetchData(&iomodel->gridonhutter,NULL,NULL,iomodel_handle,"gridonhutter"); - IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); - IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); - IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); - CreateNumberNodeToElementConnectivity(iomodel); - - for (i=0;inumberofvertices;i++){ - - /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ - if(iomodel->my_vertices[i]){ - - /*Add vertex to vertices dataset: */ - vertices->AddObject(new Vertex(i+1,i,iomodel)); - - /*Add node to nodes dataset: */ - nodes->AddObject(new Node(i+1,i,iomodel)); - - } - } - - /*Clean fetched data: */ - xfree((void**)&iomodel->deadgrids); - xfree((void**)&iomodel->x); - xfree((void**)&iomodel->y); - xfree((void**)&iomodel->z); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->gridonbed); - xfree((void**)&iomodel->gridonsurface); - xfree((void**)&iomodel->gridonhutter); - xfree((void**)&iomodel->uppernodes); - xfree((void**)&iomodel->gridonicesheet); - xfree((void**)&iomodel->gridoniceshelf); - xfree((void**)&iomodel->elements); - xfree((void**)&iomodel->numbernodetoelementconnectivity); - - /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these - * datasets, it will not be redone: */ - elements->Presort(); - nodes->Presort(); - vertices->Presort(); - materials->Presort(); - - cleanup_and_return: - - /*Assign output pointer: */ - *pelements=elements; - *pnodes=nodes; - *pvertices=vertices; - *pmaterials=materials; -} Index: /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHutter/CreateNodesDiagnosticHutter.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHutter/CreateNodesDiagnosticHutter.cpp (revision 3998) +++ /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHutter/CreateNodesDiagnosticHutter.cpp (revision 3998) @@ -0,0 +1,188 @@ +/* + * CreateElementsNodesAndMaterialsDiagnosticHutter.c: + */ + +#include "../../../DataSet/DataSet.h" +#include "../../../toolkits/toolkits.h" +#include "../../../io/io.h" +#include "../../../EnumDefinitions/EnumDefinitions.h" +#include "../../../objects/objects.h" +#include "../../../shared/shared.h" +#include "../../MeshPartitionx/MeshPartitionx.h" +#include "../../../include/include.h" +#include "../ModelProcessorx.h" + +void CreateElementsNodesAndMaterialsDiagnosticHutter(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ + + /*Intermediary*/ + int i,j,k,n; + + /*DataSets: */ + DataSet* elements = NULL; + DataSet* nodes = NULL; + DataSet* vertices = NULL; + DataSet* materials = NULL; + + /*First create the elements, nodes and material properties: */ + elements = new DataSet(ElementsEnum); + nodes = new DataSet(NodesEnum); + vertices = new DataSet(VerticesEnum); + materials = new DataSet(MaterialsEnum); + + /*Now, is the flag ishutter on? otherwise, do nothing: */ + if (!iomodel->ishutter)goto cleanup_and_return; + + /*Partition elements and vertices and nodes: */ + Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); + + /*Hutter elements can be partitioned using epart, even if + * each hutter elements either lies on a node (in 2d), or a pair of vertically juxtaposed nodes (in 3d): */ + + /*Fetch data temporarily needed: */ + IoModelFetchData(&iomodel->gridonhutter,NULL,NULL,iomodel_handle,"gridonhutter"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); + IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); + IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); + IoModelFetchData(&iomodel->drag_coefficient,NULL,NULL,iomodel_handle,"drag_coefficient"); + IoModelFetchData(&iomodel->rheology_B,NULL,NULL,iomodel_handle,"rheology_B"); + IoModelFetchData(&iomodel->rheology_n,NULL,NULL,iomodel_handle,"rheology_n"); + IoModelFetchData(&iomodel->elements_type,NULL,NULL,iomodel_handle,"elements_type"); + IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed"); + IoModelFetchData(&iomodel->elementonsurface,NULL,NULL,iomodel_handle,"elementonsurface"); + IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); + IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); + IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements"); + IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements"); + + /*2d mesh: */ + if (iomodel->dim==2){ + + for (i=0;inumberofelements;i++){ + + if(iomodel->my_elements[i]){ + + if (*(iomodel->elements_type+2*i+0)==HutterFormulationEnum){ //create only Hutter elements + + /*Create and add penta element to elements dataset: */ + elements->AddObject(new Tria(i+1,i,iomodel)); + + /*Create and add material property to materials dataset: */ + materials->AddObject(new Matice(i+1,i,iomodel,3)); + } + } + } //for (i=0;inumberofvertices;i++) + } //if (iomodel->dim==2) + else{ + + for (i=0;inumberofelements;i++){ + + if(iomodel->my_elements[i]){ + + if (*(iomodel->elements_type+2*i+0)==HutterFormulationEnum){ //create only Hutter elements + + /*Create and add penta element to elements dataset: */ + elements->AddObject(new Penta(i+1,i,iomodel)); + + + /*Create and add material property to materials dataset: */ + materials->AddObject(new Matice(i+1,i,iomodel,6)); + + } + } + + } //for (i=0;inumberofelements;i++) + + } //if (iomodel->dim==2) + + /*Free data: */ + xfree((void**)&iomodel->elements); + xfree((void**)&iomodel->elementonbed); + xfree((void**)&iomodel->elementonsurface); + xfree((void**)&iomodel->elements_type); + xfree((void**)&iomodel->gridonhutter); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->surface); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->gridonbed); + xfree((void**)&iomodel->elementonwater); + xfree((void**)&iomodel->uppernodes); + xfree((void**)&iomodel->drag_coefficient); + xfree((void**)&iomodel->rheology_B); + xfree((void**)&iomodel->rheology_n); + xfree((void**)&iomodel->upperelements); + xfree((void**)&iomodel->lowerelements); + + /*Add new constrant material property to materials, at the end: */ + if (iomodel->dim==2){ + materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel)); //put it at the end of the materials + } + else{ + materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel)); //put it at the end of the materials + } + + /*First fetch data: */ + if (iomodel->dim==3){ + IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); + IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); + } + IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); + IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); + IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); + IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); + IoModelFetchData(&iomodel->gridonhutter,NULL,NULL,iomodel_handle,"gridonhutter"); + IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); + IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); + IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); + CreateNumberNodeToElementConnectivity(iomodel); + + for (i=0;inumberofvertices;i++){ + + /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ + if(iomodel->my_vertices[i]){ + + /*Add vertex to vertices dataset: */ + vertices->AddObject(new Vertex(i+1,i,iomodel)); + + /*Add node to nodes dataset: */ + nodes->AddObject(new Node(i+1,i,iomodel)); + + } + } + + /*Clean fetched data: */ + xfree((void**)&iomodel->deadgrids); + xfree((void**)&iomodel->x); + xfree((void**)&iomodel->y); + xfree((void**)&iomodel->z); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->gridonbed); + xfree((void**)&iomodel->gridonsurface); + xfree((void**)&iomodel->gridonhutter); + xfree((void**)&iomodel->uppernodes); + xfree((void**)&iomodel->gridonicesheet); + xfree((void**)&iomodel->gridoniceshelf); + xfree((void**)&iomodel->elements); + xfree((void**)&iomodel->numbernodetoelementconnectivity); + + /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these + * datasets, it will not be redone: */ + elements->Presort(); + nodes->Presort(); + vertices->Presort(); + materials->Presort(); + + cleanup_and_return: + + /*Assign output pointer: */ + *pelements=elements; + *pnodes=nodes; + *pvertices=vertices; + *pmaterials=materials; +} Index: sm/trunk/src/c/modules/ModelProcessorx/DiagnosticStokes/CreateElementsNodesAndMaterialsDiagnosticStokes.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticStokes/CreateElementsNodesAndMaterialsDiagnosticStokes.cpp (revision 3997) +++ (revision ) @@ -1,175 +1,0 @@ -/* - * CreateElementsNodesAndMaterialsDiagnosticStokes.c: - */ - -#include "../../../DataSet/DataSet.h" -#include "../../../toolkits/toolkits.h" -#include "../../../io/io.h" -#include "../../../EnumDefinitions/EnumDefinitions.h" -#include "../../../objects/objects.h" -#include "../../../shared/shared.h" -#include "../../../include/include.h" -#include "../../MeshPartitionx/MeshPartitionx.h" -#include "../ModelProcessorx.h" - -void CreateElementsNodesAndMaterialsDiagnosticStokes(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ - - /*Intermediary*/ - int i,j,k,n; - - /*DataSets: */ - DataSet* elements = NULL; - DataSet* nodes = NULL; - DataSet* vertices = NULL; - DataSet* materials = NULL; - - /*First create the elements, nodes and material properties: */ - elements = new DataSet(ElementsEnum); - nodes = new DataSet(NodesEnum); - vertices = new DataSet(VerticesEnum); - materials = new DataSet(MaterialsEnum); - - /*Now, do we have Stokes elements?*/ - if (iomodel->dim==2) goto cleanup_and_return; - if (!iomodel->isstokes) goto cleanup_and_return; - - - /*Partition elements and vertices and nodes: */ - Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); - - /*Fetch data needed: */ - IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->drag_coefficient,NULL,NULL,iomodel_handle,"drag_coefficient"); - IoModelFetchData(&iomodel->drag_p,NULL,NULL,iomodel_handle,"drag_p"); - IoModelFetchData(&iomodel->drag_q,NULL,NULL,iomodel_handle,"drag_q"); - IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); - IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed"); - IoModelFetchData(&iomodel->elementonsurface,NULL,NULL,iomodel_handle,"elementonsurface"); - IoModelFetchData(&iomodel->elements_type,NULL,NULL,iomodel_handle,"elements_type"); - IoModelFetchData(&iomodel->rheology_B,NULL,NULL,iomodel_handle,"rheology_B"); - IoModelFetchData(&iomodel->rheology_n,NULL,NULL,iomodel_handle,"rheology_n"); - IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); - IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate"); - IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); - IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); - IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); - IoModelFetchData(&iomodel->vz,NULL,NULL,iomodel_handle,"vz"); - IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements"); - IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements"); - - if (iomodel->control_analysis){ - IoModelFetchData(&iomodel->vx_obs,NULL,NULL,iomodel_handle,"vx_obs"); - IoModelFetchData(&iomodel->vy_obs,NULL,NULL,iomodel_handle,"vy_obs"); - IoModelFetchData(&iomodel->weights,NULL,NULL,iomodel_handle,"weights"); - } - - for (i=0;inumberofelements;i++){ - - if(iomodel->my_elements[i]){ - if (*(iomodel->elements_type+2*i+1)==StokesFormulationEnum){ - - /*Create and add penta element to elements dataset: */ - elements->AddObject(new Penta(i+1,i,iomodel)); - - /*Create and add material property to materials dataset: */ - materials->AddObject(new Matice(i+1,i,iomodel,6)); - } - - } - - }//for (i=0;ielements); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->surface); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->drag_coefficient); - xfree((void**)&iomodel->drag_p); - xfree((void**)&iomodel->drag_q); - xfree((void**)&iomodel->elementoniceshelf); - xfree((void**)&iomodel->elementonbed); - xfree((void**)&iomodel->elementonsurface); - xfree((void**)&iomodel->elements_type); - xfree((void**)&iomodel->rheology_n); - xfree((void**)&iomodel->rheology_B); - xfree((void**)&iomodel->accumulation_rate); - xfree((void**)&iomodel->melting_rate); - xfree((void**)&iomodel->vx); - xfree((void**)&iomodel->vy); - xfree((void**)&iomodel->vz); - if (iomodel->control_analysis){ - xfree((void**)&iomodel->vx_obs); - xfree((void**)&iomodel->vy_obs); - xfree((void**)&iomodel->weights); - } - xfree((void**)&iomodel->upperelements); - xfree((void**)&iomodel->lowerelements); - - - /*Add new constrant material property to materials, at the end: */ - materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel)); - - /*First fetch data: */ - IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); - IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); - IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); - IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); - IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); - IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); - IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); - IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); - IoModelFetchData(&iomodel->gridonstokes,NULL,NULL,iomodel_handle,"gridonstokes"); - IoModelFetchData(&iomodel->borderstokes,NULL,NULL,iomodel_handle,"borderstokes"); - - for (i=0;inumberofvertices;i++){ - - /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ - if(iomodel->my_vertices[i]){ - - /*Add vertex to vertices dataset: */ - vertices->AddObject(new Vertex(i+1,i,iomodel)); - - /*Add node to nodes dataset: */ - nodes->AddObject(new Node(i+1,i,iomodel)); - - } - } - - /*Clean fetched data: */ - xfree((void**)&iomodel->deadgrids); - xfree((void**)&iomodel->x); - xfree((void**)&iomodel->y); - xfree((void**)&iomodel->z); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->gridonbed); - xfree((void**)&iomodel->gridonsurface); - xfree((void**)&iomodel->uppernodes); - xfree((void**)&iomodel->gridonstokes); - xfree((void**)&iomodel->borderstokes); - xfree((void**)&iomodel->gridonicesheet); - xfree((void**)&iomodel->gridoniceshelf); - - /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these - * datasets, it will not be redone: */ - elements->Presort(); - nodes->Presort(); - vertices->Presort(); - materials->Presort(); - - cleanup_and_return: - - /*Assign output pointer: */ - *pelements=elements; - *pnodes=nodes; - *pvertices=vertices; - *pmaterials=materials; - -} Index: /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticStokes/CreateNodesDiagnosticStokes.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticStokes/CreateNodesDiagnosticStokes.cpp (revision 3998) +++ /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticStokes/CreateNodesDiagnosticStokes.cpp (revision 3998) @@ -0,0 +1,175 @@ +/* + * CreateElementsNodesAndMaterialsDiagnosticStokes.c: + */ + +#include "../../../DataSet/DataSet.h" +#include "../../../toolkits/toolkits.h" +#include "../../../io/io.h" +#include "../../../EnumDefinitions/EnumDefinitions.h" +#include "../../../objects/objects.h" +#include "../../../shared/shared.h" +#include "../../../include/include.h" +#include "../../MeshPartitionx/MeshPartitionx.h" +#include "../ModelProcessorx.h" + +void CreateElementsNodesAndMaterialsDiagnosticStokes(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ + + /*Intermediary*/ + int i,j,k,n; + + /*DataSets: */ + DataSet* elements = NULL; + DataSet* nodes = NULL; + DataSet* vertices = NULL; + DataSet* materials = NULL; + + /*First create the elements, nodes and material properties: */ + elements = new DataSet(ElementsEnum); + nodes = new DataSet(NodesEnum); + vertices = new DataSet(VerticesEnum); + materials = new DataSet(MaterialsEnum); + + /*Now, do we have Stokes elements?*/ + if (iomodel->dim==2) goto cleanup_and_return; + if (!iomodel->isstokes) goto cleanup_and_return; + + + /*Partition elements and vertices and nodes: */ + Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); + + /*Fetch data needed: */ + IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->drag_coefficient,NULL,NULL,iomodel_handle,"drag_coefficient"); + IoModelFetchData(&iomodel->drag_p,NULL,NULL,iomodel_handle,"drag_p"); + IoModelFetchData(&iomodel->drag_q,NULL,NULL,iomodel_handle,"drag_q"); + IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); + IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed"); + IoModelFetchData(&iomodel->elementonsurface,NULL,NULL,iomodel_handle,"elementonsurface"); + IoModelFetchData(&iomodel->elements_type,NULL,NULL,iomodel_handle,"elements_type"); + IoModelFetchData(&iomodel->rheology_B,NULL,NULL,iomodel_handle,"rheology_B"); + IoModelFetchData(&iomodel->rheology_n,NULL,NULL,iomodel_handle,"rheology_n"); + IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); + IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate"); + IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); + IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); + IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); + IoModelFetchData(&iomodel->vz,NULL,NULL,iomodel_handle,"vz"); + IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements"); + IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements"); + + if (iomodel->control_analysis){ + IoModelFetchData(&iomodel->vx_obs,NULL,NULL,iomodel_handle,"vx_obs"); + IoModelFetchData(&iomodel->vy_obs,NULL,NULL,iomodel_handle,"vy_obs"); + IoModelFetchData(&iomodel->weights,NULL,NULL,iomodel_handle,"weights"); + } + + for (i=0;inumberofelements;i++){ + + if(iomodel->my_elements[i]){ + if (*(iomodel->elements_type+2*i+1)==StokesFormulationEnum){ + + /*Create and add penta element to elements dataset: */ + elements->AddObject(new Penta(i+1,i,iomodel)); + + /*Create and add material property to materials dataset: */ + materials->AddObject(new Matice(i+1,i,iomodel,6)); + } + + } + + }//for (i=0;ielements); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->surface); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->drag_coefficient); + xfree((void**)&iomodel->drag_p); + xfree((void**)&iomodel->drag_q); + xfree((void**)&iomodel->elementoniceshelf); + xfree((void**)&iomodel->elementonbed); + xfree((void**)&iomodel->elementonsurface); + xfree((void**)&iomodel->elements_type); + xfree((void**)&iomodel->rheology_n); + xfree((void**)&iomodel->rheology_B); + xfree((void**)&iomodel->accumulation_rate); + xfree((void**)&iomodel->melting_rate); + xfree((void**)&iomodel->vx); + xfree((void**)&iomodel->vy); + xfree((void**)&iomodel->vz); + if (iomodel->control_analysis){ + xfree((void**)&iomodel->vx_obs); + xfree((void**)&iomodel->vy_obs); + xfree((void**)&iomodel->weights); + } + xfree((void**)&iomodel->upperelements); + xfree((void**)&iomodel->lowerelements); + + + /*Add new constrant material property to materials, at the end: */ + materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel)); + + /*First fetch data: */ + IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); + IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); + IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); + IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); + IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); + IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); + IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); + IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); + IoModelFetchData(&iomodel->gridonstokes,NULL,NULL,iomodel_handle,"gridonstokes"); + IoModelFetchData(&iomodel->borderstokes,NULL,NULL,iomodel_handle,"borderstokes"); + + for (i=0;inumberofvertices;i++){ + + /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ + if(iomodel->my_vertices[i]){ + + /*Add vertex to vertices dataset: */ + vertices->AddObject(new Vertex(i+1,i,iomodel)); + + /*Add node to nodes dataset: */ + nodes->AddObject(new Node(i+1,i,iomodel)); + + } + } + + /*Clean fetched data: */ + xfree((void**)&iomodel->deadgrids); + xfree((void**)&iomodel->x); + xfree((void**)&iomodel->y); + xfree((void**)&iomodel->z); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->gridonbed); + xfree((void**)&iomodel->gridonsurface); + xfree((void**)&iomodel->uppernodes); + xfree((void**)&iomodel->gridonstokes); + xfree((void**)&iomodel->borderstokes); + xfree((void**)&iomodel->gridonicesheet); + xfree((void**)&iomodel->gridoniceshelf); + + /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these + * datasets, it will not be redone: */ + elements->Presort(); + nodes->Presort(); + vertices->Presort(); + materials->Presort(); + + cleanup_and_return: + + /*Assign output pointer: */ + *pelements=elements; + *pnodes=nodes; + *pvertices=vertices; + *pmaterials=materials; + +} Index: sm/trunk/src/c/modules/ModelProcessorx/DiagnosticVert/CreateElementsNodesAndMaterialsDiagnosticVert.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticVert/CreateElementsNodesAndMaterialsDiagnosticVert.cpp (revision 3997) +++ (revision ) @@ -1,145 +1,0 @@ -/* - * CreateElementsNodesAndMaterialsDiagnosticVert.c: - */ - -#include "../../../DataSet/DataSet.h" -#include "../../../toolkits/toolkits.h" -#include "../../../io/io.h" -#include "../../../EnumDefinitions/EnumDefinitions.h" -#include "../../../objects/objects.h" -#include "../../../shared/shared.h" -#include "../../MeshPartitionx/MeshPartitionx.h" -#include "../../../include/include.h" -#include "../ModelProcessorx.h" - -void CreateElementsNodesAndMaterialsDiagnosticVert(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ - - /*Intermediary*/ - int i,j,k,n; - - /*DataSets: */ - DataSet* elements = NULL; - DataSet* nodes = NULL; - DataSet* vertices = NULL; - DataSet* materials = NULL; - - /*First create the elements, nodes and material properties: */ - elements = new DataSet(ElementsEnum); - nodes = new DataSet(NodesEnum); - vertices = new DataSet(VerticesEnum); - materials = new DataSet(MaterialsEnum); - - /*Now, is the flag macayaealpattyn on? otherwise, do nothing: */ - if (iomodel->dim==2)goto cleanup_and_return; - - /*Partition elements and vertices and nodes: */ - Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); - - /*Create 3d elements: */ - /*Fetch data needed: */ - IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); - IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed"); - IoModelFetchData(&iomodel->elementonsurface,NULL,NULL,iomodel_handle,"elementonsurface"); - IoModelFetchData(&iomodel->elements_type,NULL,NULL,iomodel_handle,"elements_type"); - IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate"); - IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); - IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); - IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); - IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); - IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements"); - IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements"); - - for (i=0;inumberofelements;i++){ - - if(iomodel->my_elements[i]){ - - /*Create and add penta element to elements dataset: */ - elements->AddObject(new Penta(i+1,i,iomodel)); - - /*Create and add material property to materials dataset: */ - materials->AddObject(new Matice(i+1,i,iomodel,6)); - - } - - }//for (i=0;ielements); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->surface); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->elementoniceshelf); - xfree((void**)&iomodel->elementonbed); - xfree((void**)&iomodel->elementonsurface); - xfree((void**)&iomodel->elements_type); - xfree((void**)&iomodel->melting_rate); - xfree((void**)&iomodel->accumulation_rate); - xfree((void**)&iomodel->elementonwater); - xfree((void**)&iomodel->vx); - xfree((void**)&iomodel->vy); - xfree((void**)&iomodel->upperelements); - xfree((void**)&iomodel->lowerelements); - - /*Add new constrant material property to materials, at the end: */ - materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel));//put it at the end of the materials - - /*First fetch data: */ - IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); - IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); - IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); - IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); - IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); - IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); - IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); - IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); - - for (i=0;inumberofvertices;i++){ - - /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ - if(iomodel->my_vertices[i]){ - - /*Add vertex to vertices dataset: */ - vertices->AddObject(new Vertex(i+1,i,iomodel)); - - /*Add node to nodes dataset: */ - nodes->AddObject(new Node(i+1,i,iomodel)); - - } - } - - /*Clean fetched data: */ - xfree((void**)&iomodel->deadgrids); - xfree((void**)&iomodel->x); - xfree((void**)&iomodel->y); - xfree((void**)&iomodel->z); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->gridonbed); - xfree((void**)&iomodel->gridonsurface); - xfree((void**)&iomodel->uppernodes); - xfree((void**)&iomodel->gridonicesheet); - xfree((void**)&iomodel->gridoniceshelf); - - /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these - * datasets, it will not be redone: */ - elements->Presort(); - nodes->Presort(); - vertices->Presort(); - materials->Presort(); - - cleanup_and_return: - - /*Assign output pointer: */ - *pelements=elements; - *pnodes=nodes; - *pvertices=vertices; - *pmaterials=materials; - -} Index: /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticVert/CreateNodesDiagnosticVert.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticVert/CreateNodesDiagnosticVert.cpp (revision 3998) +++ /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticVert/CreateNodesDiagnosticVert.cpp (revision 3998) @@ -0,0 +1,145 @@ +/* + * CreateElementsNodesAndMaterialsDiagnosticVert.c: + */ + +#include "../../../DataSet/DataSet.h" +#include "../../../toolkits/toolkits.h" +#include "../../../io/io.h" +#include "../../../EnumDefinitions/EnumDefinitions.h" +#include "../../../objects/objects.h" +#include "../../../shared/shared.h" +#include "../../MeshPartitionx/MeshPartitionx.h" +#include "../../../include/include.h" +#include "../ModelProcessorx.h" + +void CreateElementsNodesAndMaterialsDiagnosticVert(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ + + /*Intermediary*/ + int i,j,k,n; + + /*DataSets: */ + DataSet* elements = NULL; + DataSet* nodes = NULL; + DataSet* vertices = NULL; + DataSet* materials = NULL; + + /*First create the elements, nodes and material properties: */ + elements = new DataSet(ElementsEnum); + nodes = new DataSet(NodesEnum); + vertices = new DataSet(VerticesEnum); + materials = new DataSet(MaterialsEnum); + + /*Now, is the flag macayaealpattyn on? otherwise, do nothing: */ + if (iomodel->dim==2)goto cleanup_and_return; + + /*Partition elements and vertices and nodes: */ + Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); + + /*Create 3d elements: */ + /*Fetch data needed: */ + IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); + IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed"); + IoModelFetchData(&iomodel->elementonsurface,NULL,NULL,iomodel_handle,"elementonsurface"); + IoModelFetchData(&iomodel->elements_type,NULL,NULL,iomodel_handle,"elements_type"); + IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate"); + IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); + IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); + IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); + IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); + IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements"); + IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements"); + + for (i=0;inumberofelements;i++){ + + if(iomodel->my_elements[i]){ + + /*Create and add penta element to elements dataset: */ + elements->AddObject(new Penta(i+1,i,iomodel)); + + /*Create and add material property to materials dataset: */ + materials->AddObject(new Matice(i+1,i,iomodel,6)); + + } + + }//for (i=0;ielements); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->surface); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->elementoniceshelf); + xfree((void**)&iomodel->elementonbed); + xfree((void**)&iomodel->elementonsurface); + xfree((void**)&iomodel->elements_type); + xfree((void**)&iomodel->melting_rate); + xfree((void**)&iomodel->accumulation_rate); + xfree((void**)&iomodel->elementonwater); + xfree((void**)&iomodel->vx); + xfree((void**)&iomodel->vy); + xfree((void**)&iomodel->upperelements); + xfree((void**)&iomodel->lowerelements); + + /*Add new constrant material property to materials, at the end: */ + materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel));//put it at the end of the materials + + /*First fetch data: */ + IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); + IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); + IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); + IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); + IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); + IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); + IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); + IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); + + for (i=0;inumberofvertices;i++){ + + /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ + if(iomodel->my_vertices[i]){ + + /*Add vertex to vertices dataset: */ + vertices->AddObject(new Vertex(i+1,i,iomodel)); + + /*Add node to nodes dataset: */ + nodes->AddObject(new Node(i+1,i,iomodel)); + + } + } + + /*Clean fetched data: */ + xfree((void**)&iomodel->deadgrids); + xfree((void**)&iomodel->x); + xfree((void**)&iomodel->y); + xfree((void**)&iomodel->z); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->gridonbed); + xfree((void**)&iomodel->gridonsurface); + xfree((void**)&iomodel->uppernodes); + xfree((void**)&iomodel->gridonicesheet); + xfree((void**)&iomodel->gridoniceshelf); + + /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these + * datasets, it will not be redone: */ + elements->Presort(); + nodes->Presort(); + vertices->Presort(); + materials->Presort(); + + cleanup_and_return: + + /*Assign output pointer: */ + *pelements=elements; + *pnodes=nodes; + *pvertices=vertices; + *pmaterials=materials; + +} Index: sm/trunk/src/c/modules/ModelProcessorx/Melting/CreateElementsNodesAndMaterialsMelting.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Melting/CreateElementsNodesAndMaterialsMelting.cpp (revision 3997) +++ (revision ) @@ -1,149 +1,0 @@ -/* - * CreateElementsNodesAndMaterialsMelting.c: - */ - -#include "../../../DataSet/DataSet.h" -#include "../../../toolkits/toolkits.h" -#include "../../../io/io.h" -#include "../../../EnumDefinitions/EnumDefinitions.h" -#include "../../../include/include.h" -#include "../../../objects/objects.h" -#include "../../../shared/shared.h" -#include "../../MeshPartitionx/MeshPartitionx.h" -#include "../ModelProcessorx.h" - -void CreateElementsNodesAndMaterialsMelting(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ - - /*Intermediary*/ - int i,j,k,n; - - /*DataSets: */ - DataSet* elements = NULL; - DataSet* nodes = NULL; - DataSet* vertices = NULL; - DataSet* materials = NULL; - - /*First create the elements, nodes and material properties: */ - elements = new DataSet(ElementsEnum); - nodes = new DataSet(NodesEnum); - vertices = new DataSet(VerticesEnum); - materials = new DataSet(MaterialsEnum); - - /*Partition elements and vertices and nodes: */ - Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); - - /*Fetch data needed: */ - IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->drag_coefficient,NULL,NULL,iomodel_handle,"drag_coefficient"); - IoModelFetchData(&iomodel->drag_p,NULL,NULL,iomodel_handle,"drag_p"); - IoModelFetchData(&iomodel->drag_q,NULL,NULL,iomodel_handle,"drag_q"); - IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); - IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed"); - IoModelFetchData(&iomodel->elementonsurface,NULL,NULL,iomodel_handle,"elementonsurface"); - IoModelFetchData(&iomodel->elements_type,NULL,NULL,iomodel_handle,"elements_type"); - IoModelFetchData(&iomodel->rheology_B,NULL,NULL,iomodel_handle,"rheology_B"); - IoModelFetchData(&iomodel->rheology_n,NULL,NULL,iomodel_handle,"rheology_n"); - IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); - IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate"); - IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); - IoModelFetchData(&iomodel->pressure,NULL,NULL,iomodel_handle,"pressure"); - IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements"); - IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements"); - - for (i=0;inumberofelements;i++){ - - if(iomodel->my_elements[i]){ - - /*Create and add tria element to elements dataset: */ - elements->AddObject(new Penta(i+1,i,iomodel)); - - /*Create and add material property to materials dataset: */ - materials->AddObject(new Matice(i+1,i,iomodel,6)); - } - - }//for (i=0;ielements); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->surface); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->drag_coefficient); - xfree((void**)&iomodel->drag_p); - xfree((void**)&iomodel->drag_q); - xfree((void**)&iomodel->elementoniceshelf); - xfree((void**)&iomodel->elementonbed); - xfree((void**)&iomodel->elementonsurface); - xfree((void**)&iomodel->elements_type); - xfree((void**)&iomodel->rheology_n); - xfree((void**)&iomodel->rheology_B); - xfree((void**)&iomodel->accumulation_rate); - xfree((void**)&iomodel->melting_rate); - xfree((void**)&iomodel->elementonwater); - xfree((void**)&iomodel->pressure); - xfree((void**)&iomodel->upperelements); - xfree((void**)&iomodel->lowerelements); - - /*Add new constrant material property tgo materials, at the end: */ - materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel));//put it at the end of the materials - - /*First fetch data: */ - if (iomodel->dim==3){ - IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); - IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); - } - IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); - IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); - IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); - IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); - IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); - IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); - - for (i=0;inumberofvertices;i++){ - - /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ - if(iomodel->my_vertices[i]){ - - /*Add vertex to vertices dataset: */ - vertices->AddObject(new Vertex(i+1,i,iomodel)); - - /*Add node to nodes dataset: */ - nodes->AddObject(new Node(i+1,i,iomodel)); - - } - } - - /*Clean fetched data: */ - xfree((void**)&iomodel->deadgrids); - xfree((void**)&iomodel->x); - xfree((void**)&iomodel->y); - xfree((void**)&iomodel->z); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->gridonbed); - xfree((void**)&iomodel->gridonsurface); - xfree((void**)&iomodel->uppernodes); - xfree((void**)&iomodel->gridonicesheet); - xfree((void**)&iomodel->gridoniceshelf); - - /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these - * datasets, it will not be redone: */ - elements->Presort(); - nodes->Presort(); - vertices->Presort(); - materials->Presort(); - - cleanup_and_return: - - /*Assign output pointer: */ - *pelements=elements; - *pnodes=nodes; - *pvertices=vertices; - *pmaterials=materials; -} Index: /issm/trunk/src/c/modules/ModelProcessorx/Melting/CreateNodesMelting.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Melting/CreateNodesMelting.cpp (revision 3998) +++ /issm/trunk/src/c/modules/ModelProcessorx/Melting/CreateNodesMelting.cpp (revision 3998) @@ -0,0 +1,149 @@ +/* + * CreateElementsNodesAndMaterialsMelting.c: + */ + +#include "../../../DataSet/DataSet.h" +#include "../../../toolkits/toolkits.h" +#include "../../../io/io.h" +#include "../../../EnumDefinitions/EnumDefinitions.h" +#include "../../../include/include.h" +#include "../../../objects/objects.h" +#include "../../../shared/shared.h" +#include "../../MeshPartitionx/MeshPartitionx.h" +#include "../ModelProcessorx.h" + +void CreateElementsNodesAndMaterialsMelting(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ + + /*Intermediary*/ + int i,j,k,n; + + /*DataSets: */ + DataSet* elements = NULL; + DataSet* nodes = NULL; + DataSet* vertices = NULL; + DataSet* materials = NULL; + + /*First create the elements, nodes and material properties: */ + elements = new DataSet(ElementsEnum); + nodes = new DataSet(NodesEnum); + vertices = new DataSet(VerticesEnum); + materials = new DataSet(MaterialsEnum); + + /*Partition elements and vertices and nodes: */ + Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); + + /*Fetch data needed: */ + IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->drag_coefficient,NULL,NULL,iomodel_handle,"drag_coefficient"); + IoModelFetchData(&iomodel->drag_p,NULL,NULL,iomodel_handle,"drag_p"); + IoModelFetchData(&iomodel->drag_q,NULL,NULL,iomodel_handle,"drag_q"); + IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); + IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed"); + IoModelFetchData(&iomodel->elementonsurface,NULL,NULL,iomodel_handle,"elementonsurface"); + IoModelFetchData(&iomodel->elements_type,NULL,NULL,iomodel_handle,"elements_type"); + IoModelFetchData(&iomodel->rheology_B,NULL,NULL,iomodel_handle,"rheology_B"); + IoModelFetchData(&iomodel->rheology_n,NULL,NULL,iomodel_handle,"rheology_n"); + IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); + IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate"); + IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); + IoModelFetchData(&iomodel->pressure,NULL,NULL,iomodel_handle,"pressure"); + IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements"); + IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements"); + + for (i=0;inumberofelements;i++){ + + if(iomodel->my_elements[i]){ + + /*Create and add tria element to elements dataset: */ + elements->AddObject(new Penta(i+1,i,iomodel)); + + /*Create and add material property to materials dataset: */ + materials->AddObject(new Matice(i+1,i,iomodel,6)); + } + + }//for (i=0;ielements); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->surface); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->drag_coefficient); + xfree((void**)&iomodel->drag_p); + xfree((void**)&iomodel->drag_q); + xfree((void**)&iomodel->elementoniceshelf); + xfree((void**)&iomodel->elementonbed); + xfree((void**)&iomodel->elementonsurface); + xfree((void**)&iomodel->elements_type); + xfree((void**)&iomodel->rheology_n); + xfree((void**)&iomodel->rheology_B); + xfree((void**)&iomodel->accumulation_rate); + xfree((void**)&iomodel->melting_rate); + xfree((void**)&iomodel->elementonwater); + xfree((void**)&iomodel->pressure); + xfree((void**)&iomodel->upperelements); + xfree((void**)&iomodel->lowerelements); + + /*Add new constrant material property tgo materials, at the end: */ + materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel));//put it at the end of the materials + + /*First fetch data: */ + if (iomodel->dim==3){ + IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); + IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); + } + IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); + IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); + IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); + IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); + IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); + IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); + + for (i=0;inumberofvertices;i++){ + + /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ + if(iomodel->my_vertices[i]){ + + /*Add vertex to vertices dataset: */ + vertices->AddObject(new Vertex(i+1,i,iomodel)); + + /*Add node to nodes dataset: */ + nodes->AddObject(new Node(i+1,i,iomodel)); + + } + } + + /*Clean fetched data: */ + xfree((void**)&iomodel->deadgrids); + xfree((void**)&iomodel->x); + xfree((void**)&iomodel->y); + xfree((void**)&iomodel->z); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->gridonbed); + xfree((void**)&iomodel->gridonsurface); + xfree((void**)&iomodel->uppernodes); + xfree((void**)&iomodel->gridonicesheet); + xfree((void**)&iomodel->gridoniceshelf); + + /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these + * datasets, it will not be redone: */ + elements->Presort(); + nodes->Presort(); + vertices->Presort(); + materials->Presort(); + + cleanup_and_return: + + /*Assign output pointer: */ + *pelements=elements; + *pnodes=nodes; + *pvertices=vertices; + *pmaterials=materials; +} Index: sm/trunk/src/c/modules/ModelProcessorx/Prognostic/CreateElementsNodesAndMaterialsPrognostic.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Prognostic/CreateElementsNodesAndMaterialsPrognostic.cpp (revision 3997) +++ (revision ) @@ -1,187 +1,0 @@ -/* - * CreateElementsNodesAndMaterialsPrognostic.c: - */ - -#include "../../../DataSet/DataSet.h" -#include "../../../toolkits/toolkits.h" -#include "../../../EnumDefinitions/EnumDefinitions.h" -#include "../../../objects/objects.h" -#include "../../../shared/shared.h" -#include "../../../io/io.h" -#include "../../MeshPartitionx/MeshPartitionx.h" -#include "../../../include/include.h" -#include "../ModelProcessorx.h" - -void CreateElementsNodesAndMaterialsPrognostic(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ - - /*Intermediary*/ - int i,j,k,n; - - /*DataSets: */ - DataSet* elements = NULL; - DataSet* nodes = NULL; - DataSet* vertices = NULL; - DataSet* materials = NULL; - - /*First create the elements, nodes and material properties: */ - elements = new DataSet(ElementsEnum); - nodes = new DataSet(NodesEnum); - vertices = new DataSet(VerticesEnum); - materials = new DataSet(MaterialsEnum); - - /*Partition elements and vertices and nodes: */ - Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); - - /*2d mesh: */ - if (iomodel->dim==2){ - - /*Fetch data needed: */ - IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); - IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); - IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate"); - IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); - IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); - IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); - - for (i=0;inumberofelements;i++){ - - if(iomodel->my_elements[i]){ - - /*Create and add tria element to elements dataset: */ - elements->AddObject(new Tria(i+1,i,iomodel)); - - /*Create and add material property to materials dataset: */ - materials->AddObject(new Matice(i+1,i,iomodel,3)); - } - }//for (i=0;ielements); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->surface); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->elementoniceshelf); - xfree((void**)&iomodel->elementonwater); - xfree((void**)&iomodel->melting_rate); - xfree((void**)&iomodel->accumulation_rate); - xfree((void**)&iomodel->vx); - xfree((void**)&iomodel->vy); - - } - else{ // if (dim==2) - - /*Fetch data needed: */ - IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); - IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed"); - IoModelFetchData(&iomodel->elementonsurface,NULL,NULL,iomodel_handle,"elementonsurface"); - IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); - IoModelFetchData(&iomodel->pressure,NULL,NULL,iomodel_handle,"pressure"); - IoModelFetchData(&iomodel->temperature,NULL,NULL,iomodel_handle,"temperature"); - IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate"); - IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); - IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); - IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); - IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements"); - IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements"); - - for (i=0;inumberofelements;i++){ - if(iomodel->my_elements[i]){ - /*Create and add penta element to elements dataset: */ - elements->AddObject(new Penta(i+1,i,iomodel)); - - /*Create and add material property to materials dataset: */ - materials->AddObject(new Matice(i+1,i,iomodel,6)); - } - }//for (i=0;ielements); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->surface); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->elementoniceshelf); - xfree((void**)&iomodel->elementonbed); - xfree((void**)&iomodel->elementonsurface); - xfree((void**)&iomodel->elementonwater); - xfree((void**)&iomodel->pressure); - xfree((void**)&iomodel->temperature); - xfree((void**)&iomodel->melting_rate); - xfree((void**)&iomodel->accumulation_rate); - xfree((void**)&iomodel->vx); - xfree((void**)&iomodel->vy); - xfree((void**)&iomodel->upperelements); - xfree((void**)&iomodel->lowerelements); - - - } //if (dim==2) - - /*Add new constrant material property to materials, at the end: */ - materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel));//put it at the end of the materials - - /*First fetch data: */ - if (iomodel->dim==3){ - IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); - IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); - } - IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); - IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); - IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); - IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); - IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); - IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); - - for (i=0;inumberofvertices;i++){ - - /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ - if(iomodel->my_vertices[i]){ - - /*Add vertex to vertices dataset: */ - vertices->AddObject(new Vertex(i+1,i,iomodel)); - - /*Add node to nodes dataset: */ - nodes->AddObject(new Node(i+1,i,iomodel)); - - } - } - - /*Clean fetched data: */ - xfree((void**)&iomodel->deadgrids); - xfree((void**)&iomodel->x); - xfree((void**)&iomodel->y); - xfree((void**)&iomodel->z); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->gridonbed); - xfree((void**)&iomodel->gridonsurface); - xfree((void**)&iomodel->uppernodes); - xfree((void**)&iomodel->gridonicesheet); - xfree((void**)&iomodel->gridoniceshelf); - - /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these - * datasets, it will not be redone: */ - elements->Presort(); - nodes->Presort(); - vertices->Presort(); - materials->Presort(); - -cleanup_and_return: - - /*Assign output pointer: */ - *pelements=elements; - *pnodes=nodes; - *pvertices=vertices; - *pmaterials=materials; - -} Index: /issm/trunk/src/c/modules/ModelProcessorx/Prognostic/CreateNodesPrognostic.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Prognostic/CreateNodesPrognostic.cpp (revision 3998) +++ /issm/trunk/src/c/modules/ModelProcessorx/Prognostic/CreateNodesPrognostic.cpp (revision 3998) @@ -0,0 +1,187 @@ +/* + * CreateElementsNodesAndMaterialsPrognostic.c: + */ + +#include "../../../DataSet/DataSet.h" +#include "../../../toolkits/toolkits.h" +#include "../../../EnumDefinitions/EnumDefinitions.h" +#include "../../../objects/objects.h" +#include "../../../shared/shared.h" +#include "../../../io/io.h" +#include "../../MeshPartitionx/MeshPartitionx.h" +#include "../../../include/include.h" +#include "../ModelProcessorx.h" + +void CreateElementsNodesAndMaterialsPrognostic(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ + + /*Intermediary*/ + int i,j,k,n; + + /*DataSets: */ + DataSet* elements = NULL; + DataSet* nodes = NULL; + DataSet* vertices = NULL; + DataSet* materials = NULL; + + /*First create the elements, nodes and material properties: */ + elements = new DataSet(ElementsEnum); + nodes = new DataSet(NodesEnum); + vertices = new DataSet(VerticesEnum); + materials = new DataSet(MaterialsEnum); + + /*Partition elements and vertices and nodes: */ + Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); + + /*2d mesh: */ + if (iomodel->dim==2){ + + /*Fetch data needed: */ + IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); + IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); + IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate"); + IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); + IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); + IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); + + for (i=0;inumberofelements;i++){ + + if(iomodel->my_elements[i]){ + + /*Create and add tria element to elements dataset: */ + elements->AddObject(new Tria(i+1,i,iomodel)); + + /*Create and add material property to materials dataset: */ + materials->AddObject(new Matice(i+1,i,iomodel,3)); + } + }//for (i=0;ielements); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->surface); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->elementoniceshelf); + xfree((void**)&iomodel->elementonwater); + xfree((void**)&iomodel->melting_rate); + xfree((void**)&iomodel->accumulation_rate); + xfree((void**)&iomodel->vx); + xfree((void**)&iomodel->vy); + + } + else{ // if (dim==2) + + /*Fetch data needed: */ + IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); + IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed"); + IoModelFetchData(&iomodel->elementonsurface,NULL,NULL,iomodel_handle,"elementonsurface"); + IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); + IoModelFetchData(&iomodel->pressure,NULL,NULL,iomodel_handle,"pressure"); + IoModelFetchData(&iomodel->temperature,NULL,NULL,iomodel_handle,"temperature"); + IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate"); + IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); + IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); + IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); + IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements"); + IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements"); + + for (i=0;inumberofelements;i++){ + if(iomodel->my_elements[i]){ + /*Create and add penta element to elements dataset: */ + elements->AddObject(new Penta(i+1,i,iomodel)); + + /*Create and add material property to materials dataset: */ + materials->AddObject(new Matice(i+1,i,iomodel,6)); + } + }//for (i=0;ielements); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->surface); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->elementoniceshelf); + xfree((void**)&iomodel->elementonbed); + xfree((void**)&iomodel->elementonsurface); + xfree((void**)&iomodel->elementonwater); + xfree((void**)&iomodel->pressure); + xfree((void**)&iomodel->temperature); + xfree((void**)&iomodel->melting_rate); + xfree((void**)&iomodel->accumulation_rate); + xfree((void**)&iomodel->vx); + xfree((void**)&iomodel->vy); + xfree((void**)&iomodel->upperelements); + xfree((void**)&iomodel->lowerelements); + + + } //if (dim==2) + + /*Add new constrant material property to materials, at the end: */ + materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel));//put it at the end of the materials + + /*First fetch data: */ + if (iomodel->dim==3){ + IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); + IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); + } + IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); + IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); + IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); + IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); + IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); + IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); + + for (i=0;inumberofvertices;i++){ + + /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ + if(iomodel->my_vertices[i]){ + + /*Add vertex to vertices dataset: */ + vertices->AddObject(new Vertex(i+1,i,iomodel)); + + /*Add node to nodes dataset: */ + nodes->AddObject(new Node(i+1,i,iomodel)); + + } + } + + /*Clean fetched data: */ + xfree((void**)&iomodel->deadgrids); + xfree((void**)&iomodel->x); + xfree((void**)&iomodel->y); + xfree((void**)&iomodel->z); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->gridonbed); + xfree((void**)&iomodel->gridonsurface); + xfree((void**)&iomodel->uppernodes); + xfree((void**)&iomodel->gridonicesheet); + xfree((void**)&iomodel->gridoniceshelf); + + /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these + * datasets, it will not be redone: */ + elements->Presort(); + nodes->Presort(); + vertices->Presort(); + materials->Presort(); + +cleanup_and_return: + + /*Assign output pointer: */ + *pelements=elements; + *pnodes=nodes; + *pvertices=vertices; + *pmaterials=materials; + +} Index: sm/trunk/src/c/modules/ModelProcessorx/Prognostic2/CreateElementsNodesAndMaterialsPrognostic2.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Prognostic2/CreateElementsNodesAndMaterialsPrognostic2.cpp (revision 3997) +++ (revision ) @@ -1,160 +1,0 @@ -/* - * CreateElementsNodesAndMaterialsPrognostic2.c: - */ - -#include "../../../DataSet/DataSet.h" -#include "../../../toolkits/toolkits.h" -#include "../../../EnumDefinitions/EnumDefinitions.h" -#include "../../../objects/objects.h" -#include "../../../io/io.h" -#include "../../../shared/shared.h" -#include "../../../include/include.h" -#include "../../MeshPartitionx/MeshPartitionx.h" -#include "../ModelProcessorx.h" - -void CreateElementsNodesAndMaterialsPrognostic2(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ - - /*Intermediary*/ - int i,j; - int vertex_index; - int node_index; - - /*DataSets: */ - DataSet* elements = NULL; - DataSet* nodes = NULL; - DataSet* vertices = NULL; - DataSet* materials = NULL; - - /*First create the elements, nodes and material properties: */ - elements = new DataSet(ElementsEnum); - nodes = new DataSet(NodesEnum); - vertices = new DataSet(VerticesEnum); - materials = new DataSet(MaterialsEnum); - - /*Partition elements and vertices and nodes: */ - Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); - - /*elements created vary if we are dealing with a 2d mesh, or a 3d mesh: */ - /*2d mesh: */ - if (iomodel->dim==2){ - - /*Fetch data needed: */ - IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); - IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); - IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate"); - IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); - IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); - IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); - - for (i=0;inumberofelements;i++){ - - if(iomodel->my_elements[i]){ - - /*Create and add tria element to elements dataset: */ - elements->AddObject(new Tria(i+1,i,iomodel)); - - /*Create and add material property to materials dataset: */ - materials->AddObject(new Matice(i+1,i,iomodel,3)); - } - }//for (i=0;ithickness); - xfree((void**)&iomodel->surface); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->elementoniceshelf); - xfree((void**)&iomodel->elementonwater); - xfree((void**)&iomodel->melting_rate); - xfree((void**)&iomodel->accumulation_rate); - xfree((void**)&iomodel->vx); - xfree((void**)&iomodel->vy); - - } - else{ // if (dim==2) - ISSMERROR("not implemented yet"); - } //if (dim==2) - - /*Add new constrant material property tgo materials, at the end: */ - materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel)); - - /*Create nodes and vertices: */ - IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); - IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); - IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); - IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); - IoModelFetchData(&iomodel->gridonhutter,NULL,NULL,iomodel_handle,"gridonhutter"); - IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); - IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); - if (iomodel->dim==3){ - IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); - IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); - } - - /*Build Vertices dataset*/ - for (i=0;inumberofvertices;i++){ - - /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ - if(iomodel->my_vertices[i]){ - - /*Add vertex to vertices dataset: */ - vertices->AddObject(new Vertex(i+1,i,iomodel)); - - } - } - - /*Build Nodes dataset -> 3 for each element (Discontinuous Galerkin)*/ - for (i=0;inumberofelements;i++){ - for (j=0;j<3;j++){ - - if(iomodel->my_nodes[3*i+j]){ - - //Get index of the vertex on which the current node is located - vertex_index=(int)*(iomodel->elements+3*i+j)-1; //(Matlab to C indexing) - ISSMASSERT(vertex_index>=0 && vertex_indexnumberofvertices); - - //Compute Node index (id-1) - node_index=3*i+j; - - /*Add node to nodes dataset: */ - nodes->AddObject(new Node(node_index+1,vertex_index,node_index,iomodel)); - - } - } - } - - /*Clean fetched data: */ - xfree((void**)&iomodel->deadgrids); - xfree((void**)&iomodel->x); - xfree((void**)&iomodel->y); - xfree((void**)&iomodel->z); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->gridonbed); - xfree((void**)&iomodel->gridonsurface); - xfree((void**)&iomodel->gridonhutter); - xfree((void**)&iomodel->uppernodes); - xfree((void**)&iomodel->gridonicesheet); - xfree((void**)&iomodel->gridoniceshelf); - - /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these - * datasets, it will not be redone: */ - elements->Presort(); - nodes->Presort(); - vertices->Presort(); - materials->Presort(); - - cleanup_and_return: - - /*Assign output pointer: */ - *pelements=elements; - *pnodes=nodes; - *pvertices=vertices; - *pmaterials=materials; -} Index: /issm/trunk/src/c/modules/ModelProcessorx/Prognostic2/CreateNodesPrognostic2.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Prognostic2/CreateNodesPrognostic2.cpp (revision 3998) +++ /issm/trunk/src/c/modules/ModelProcessorx/Prognostic2/CreateNodesPrognostic2.cpp (revision 3998) @@ -0,0 +1,160 @@ +/* + * CreateElementsNodesAndMaterialsPrognostic2.c: + */ + +#include "../../../DataSet/DataSet.h" +#include "../../../toolkits/toolkits.h" +#include "../../../EnumDefinitions/EnumDefinitions.h" +#include "../../../objects/objects.h" +#include "../../../io/io.h" +#include "../../../shared/shared.h" +#include "../../../include/include.h" +#include "../../MeshPartitionx/MeshPartitionx.h" +#include "../ModelProcessorx.h" + +void CreateElementsNodesAndMaterialsPrognostic2(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ + + /*Intermediary*/ + int i,j; + int vertex_index; + int node_index; + + /*DataSets: */ + DataSet* elements = NULL; + DataSet* nodes = NULL; + DataSet* vertices = NULL; + DataSet* materials = NULL; + + /*First create the elements, nodes and material properties: */ + elements = new DataSet(ElementsEnum); + nodes = new DataSet(NodesEnum); + vertices = new DataSet(VerticesEnum); + materials = new DataSet(MaterialsEnum); + + /*Partition elements and vertices and nodes: */ + Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); + + /*elements created vary if we are dealing with a 2d mesh, or a 3d mesh: */ + /*2d mesh: */ + if (iomodel->dim==2){ + + /*Fetch data needed: */ + IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); + IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); + IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate"); + IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate"); + IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); + IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); + + for (i=0;inumberofelements;i++){ + + if(iomodel->my_elements[i]){ + + /*Create and add tria element to elements dataset: */ + elements->AddObject(new Tria(i+1,i,iomodel)); + + /*Create and add material property to materials dataset: */ + materials->AddObject(new Matice(i+1,i,iomodel,3)); + } + }//for (i=0;ithickness); + xfree((void**)&iomodel->surface); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->elementoniceshelf); + xfree((void**)&iomodel->elementonwater); + xfree((void**)&iomodel->melting_rate); + xfree((void**)&iomodel->accumulation_rate); + xfree((void**)&iomodel->vx); + xfree((void**)&iomodel->vy); + + } + else{ // if (dim==2) + ISSMERROR("not implemented yet"); + } //if (dim==2) + + /*Add new constrant material property tgo materials, at the end: */ + materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel)); + + /*Create nodes and vertices: */ + IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); + IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); + IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); + IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); + IoModelFetchData(&iomodel->gridonhutter,NULL,NULL,iomodel_handle,"gridonhutter"); + IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); + IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); + if (iomodel->dim==3){ + IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); + IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); + } + + /*Build Vertices dataset*/ + for (i=0;inumberofvertices;i++){ + + /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ + if(iomodel->my_vertices[i]){ + + /*Add vertex to vertices dataset: */ + vertices->AddObject(new Vertex(i+1,i,iomodel)); + + } + } + + /*Build Nodes dataset -> 3 for each element (Discontinuous Galerkin)*/ + for (i=0;inumberofelements;i++){ + for (j=0;j<3;j++){ + + if(iomodel->my_nodes[3*i+j]){ + + //Get index of the vertex on which the current node is located + vertex_index=(int)*(iomodel->elements+3*i+j)-1; //(Matlab to C indexing) + ISSMASSERT(vertex_index>=0 && vertex_indexnumberofvertices); + + //Compute Node index (id-1) + node_index=3*i+j; + + /*Add node to nodes dataset: */ + nodes->AddObject(new Node(node_index+1,vertex_index,node_index,iomodel)); + + } + } + } + + /*Clean fetched data: */ + xfree((void**)&iomodel->deadgrids); + xfree((void**)&iomodel->x); + xfree((void**)&iomodel->y); + xfree((void**)&iomodel->z); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->gridonbed); + xfree((void**)&iomodel->gridonsurface); + xfree((void**)&iomodel->gridonhutter); + xfree((void**)&iomodel->uppernodes); + xfree((void**)&iomodel->gridonicesheet); + xfree((void**)&iomodel->gridoniceshelf); + + /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these + * datasets, it will not be redone: */ + elements->Presort(); + nodes->Presort(); + vertices->Presort(); + materials->Presort(); + + cleanup_and_return: + + /*Assign output pointer: */ + *pelements=elements; + *pnodes=nodes; + *pvertices=vertices; + *pmaterials=materials; +} Index: sm/trunk/src/c/modules/ModelProcessorx/SlopeCompute/CreateElementsNodesAndMaterialsSlopeCompute.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/SlopeCompute/CreateElementsNodesAndMaterialsSlopeCompute.cpp (revision 3997) +++ (revision ) @@ -1,158 +1,0 @@ -/* - * CreateElementsNodesAndMaterialsSlopeCompute.c: - */ - -#include "../../../DataSet/DataSet.h" -#include "../../../io/io.h" -#include "../../../toolkits/toolkits.h" -#include "../../../EnumDefinitions/EnumDefinitions.h" -#include "../../../objects/objects.h" -#include "../../../shared/shared.h" -#include "../../../include/include.h" -#include "../../MeshPartitionx/MeshPartitionx.h" -#include "../ModelProcessorx.h" - -void CreateElementsNodesAndMaterialsSlopeCompute(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ - - /*output: int* epart, int* my_grids, double* my_bordergrids*/ - - int i,j,k,n; - /*DataSets: */ - DataSet* elements = NULL; - DataSet* nodes = NULL; - DataSet* vertices = NULL; - DataSet* materials = NULL; - - /*First create the elements, nodes and material properties: */ - elements = new DataSet(ElementsEnum); - nodes = new DataSet(NodesEnum); - vertices = new DataSet(VerticesEnum); - materials = new DataSet(MaterialsEnum); - - /*Partition elements and vertices and nodes: */ - Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); - - /*2d mesh: */ - if (iomodel->dim==2){ - - /*Fetch data needed: */ - IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); - IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); - IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - - for (i=0;inumberofelements;i++){ - - if(iomodel->my_elements[i]){ - - /*Create and add tria element to elements dataset: */ - elements->AddObject(new Tria(i+1,i,iomodel)); - - /*Create and add material property to materials dataset: */ - materials->AddObject(new Matice(i+1,i,iomodel,3)); - } - - }//for (i=0;ielements); - xfree((void**)&iomodel->surface); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->elementonwater); - - } - else{ // if (dim==2) - - /*Fetch data needed: */ - IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); - IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed"); - IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); - IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements"); - IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements"); - - for (i=0;inumberofelements;i++){ - if(iomodel->my_elements[i]){ - /*Create and add penta element to elements dataset: */ - elements->AddObject(new Penta(i+1,i,iomodel)); - - /*Create and add material property to materials dataset: */ - materials->AddObject(new Matice(i+1,i,iomodel,6)); - - }//if(my_elements[i]) - }//for (i=0;ielements); - xfree((void**)&iomodel->surface); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->elementonbed); - xfree((void**)&iomodel->elementonwater); - xfree((void**)&iomodel->upperelements); - xfree((void**)&iomodel->lowerelements); - - } //if (dim==2) - - /*Add new constrant material property tgo materials, at the end: */ - materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel)); //put it at the end of the materials - - /*First fetch data: */ - if (iomodel->dim==3){ - IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); - IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); - } - IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); - IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); - IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); - IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); - IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); - IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); - - for (i=0;inumberofvertices;i++){ - - /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ - if(iomodel->my_vertices[i]){ - - /*Add vertex to vertices dataset: */ - vertices->AddObject(new Vertex(i+1,i,iomodel)); - - /*Add node to nodes dataset: */ - nodes->AddObject(new Node(i+1,i,iomodel)); - - } - } - - /*Clean fetched data: */ - xfree((void**)&iomodel->deadgrids); - xfree((void**)&iomodel->x); - xfree((void**)&iomodel->y); - xfree((void**)&iomodel->z); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->gridonbed); - xfree((void**)&iomodel->gridonsurface); - xfree((void**)&iomodel->uppernodes); - xfree((void**)&iomodel->gridonicesheet); - xfree((void**)&iomodel->gridoniceshelf); - - /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these - * datasets, it will not be redone: */ - elements->Presort(); - nodes->Presort(); - vertices->Presort(); - materials->Presort(); - - cleanup_and_return: - - /*Assign output pointer: */ - *pelements=elements; - *pnodes=nodes; - *pvertices=vertices; - *pmaterials=materials; - -} Index: /issm/trunk/src/c/modules/ModelProcessorx/SlopeCompute/CreateNodesSlopeConpute.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/SlopeCompute/CreateNodesSlopeConpute.cpp (revision 3998) +++ /issm/trunk/src/c/modules/ModelProcessorx/SlopeCompute/CreateNodesSlopeConpute.cpp (revision 3998) @@ -0,0 +1,158 @@ +/* + * CreateElementsNodesAndMaterialsSlopeCompute.c: + */ + +#include "../../../DataSet/DataSet.h" +#include "../../../io/io.h" +#include "../../../toolkits/toolkits.h" +#include "../../../EnumDefinitions/EnumDefinitions.h" +#include "../../../objects/objects.h" +#include "../../../shared/shared.h" +#include "../../../include/include.h" +#include "../../MeshPartitionx/MeshPartitionx.h" +#include "../ModelProcessorx.h" + +void CreateElementsNodesAndMaterialsSlopeCompute(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ + + /*output: int* epart, int* my_grids, double* my_bordergrids*/ + + int i,j,k,n; + /*DataSets: */ + DataSet* elements = NULL; + DataSet* nodes = NULL; + DataSet* vertices = NULL; + DataSet* materials = NULL; + + /*First create the elements, nodes and material properties: */ + elements = new DataSet(ElementsEnum); + nodes = new DataSet(NodesEnum); + vertices = new DataSet(VerticesEnum); + materials = new DataSet(MaterialsEnum); + + /*Partition elements and vertices and nodes: */ + Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); + + /*2d mesh: */ + if (iomodel->dim==2){ + + /*Fetch data needed: */ + IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); + IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); + IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + + for (i=0;inumberofelements;i++){ + + if(iomodel->my_elements[i]){ + + /*Create and add tria element to elements dataset: */ + elements->AddObject(new Tria(i+1,i,iomodel)); + + /*Create and add material property to materials dataset: */ + materials->AddObject(new Matice(i+1,i,iomodel,3)); + } + + }//for (i=0;ielements); + xfree((void**)&iomodel->surface); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->elementonwater); + + } + else{ // if (dim==2) + + /*Fetch data needed: */ + IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); + IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed"); + IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); + IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements"); + IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements"); + + for (i=0;inumberofelements;i++){ + if(iomodel->my_elements[i]){ + /*Create and add penta element to elements dataset: */ + elements->AddObject(new Penta(i+1,i,iomodel)); + + /*Create and add material property to materials dataset: */ + materials->AddObject(new Matice(i+1,i,iomodel,6)); + + }//if(my_elements[i]) + }//for (i=0;ielements); + xfree((void**)&iomodel->surface); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->elementonbed); + xfree((void**)&iomodel->elementonwater); + xfree((void**)&iomodel->upperelements); + xfree((void**)&iomodel->lowerelements); + + } //if (dim==2) + + /*Add new constrant material property tgo materials, at the end: */ + materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel)); //put it at the end of the materials + + /*First fetch data: */ + if (iomodel->dim==3){ + IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); + IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); + } + IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); + IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); + IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); + IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); + IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); + IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); + + for (i=0;inumberofvertices;i++){ + + /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ + if(iomodel->my_vertices[i]){ + + /*Add vertex to vertices dataset: */ + vertices->AddObject(new Vertex(i+1,i,iomodel)); + + /*Add node to nodes dataset: */ + nodes->AddObject(new Node(i+1,i,iomodel)); + + } + } + + /*Clean fetched data: */ + xfree((void**)&iomodel->deadgrids); + xfree((void**)&iomodel->x); + xfree((void**)&iomodel->y); + xfree((void**)&iomodel->z); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->gridonbed); + xfree((void**)&iomodel->gridonsurface); + xfree((void**)&iomodel->uppernodes); + xfree((void**)&iomodel->gridonicesheet); + xfree((void**)&iomodel->gridoniceshelf); + + /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these + * datasets, it will not be redone: */ + elements->Presort(); + nodes->Presort(); + vertices->Presort(); + materials->Presort(); + + cleanup_and_return: + + /*Assign output pointer: */ + *pelements=elements; + *pnodes=nodes; + *pvertices=vertices; + *pmaterials=materials; + +} Index: sm/trunk/src/c/modules/ModelProcessorx/Thermal/CreateElementsNodesAndMaterialsThermal.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Thermal/CreateElementsNodesAndMaterialsThermal.cpp (revision 3997) +++ (revision ) @@ -1,154 +1,0 @@ -/* - * CreateElementsNodesAndMaterialsThermal.c: - */ - -#include "../../../DataSet/DataSet.h" -#include "../../../toolkits/toolkits.h" -#include "../../../io/io.h" -#include "../../../EnumDefinitions/EnumDefinitions.h" -#include "../../../objects/objects.h" -#include "../../../include/include.h" -#include "../../../shared/shared.h" -#include "../../MeshPartitionx/MeshPartitionx.h" -#include "../ModelProcessorx.h" - -void CreateElementsNodesAndMaterialsThermal(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ - - /*Intermediary*/ - int i,j,k,n; - - /*DataSets: */ - DataSet* elements = NULL; - DataSet* nodes = NULL; - DataSet* vertices = NULL; - DataSet* materials = NULL; - - /*First create the elements, nodes and material properties: */ - elements = new DataSet(ElementsEnum); - nodes = new DataSet(NodesEnum); - vertices = new DataSet(VerticesEnum); - materials = new DataSet(MaterialsEnum); - - /*Partition elements and vertices and nodes: */ - Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); - - /*Fetch data needed: */ - IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->drag_coefficient,NULL,NULL,iomodel_handle,"drag_coefficient"); - IoModelFetchData(&iomodel->drag_p,NULL,NULL,iomodel_handle,"drag_p"); - IoModelFetchData(&iomodel->drag_q,NULL,NULL,iomodel_handle,"drag_q"); - IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); - IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed"); - IoModelFetchData(&iomodel->elementonsurface,NULL,NULL,iomodel_handle,"elementonsurface"); - IoModelFetchData(&iomodel->elements_type,NULL,NULL,iomodel_handle,"elements_type"); - IoModelFetchData(&iomodel->geothermalflux,NULL,NULL,iomodel_handle,"geothermalflux"); - IoModelFetchData(&iomodel->rheology_B,NULL,NULL,iomodel_handle,"rheology_B"); - IoModelFetchData(&iomodel->rheology_n,NULL,NULL,iomodel_handle,"rheology_n"); - IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); - IoModelFetchData(&iomodel->pressure,NULL,NULL,iomodel_handle,"pressure"); - IoModelFetchData(&iomodel->temperature,NULL,NULL,iomodel_handle,"temperature"); - IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); - IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); - IoModelFetchData(&iomodel->vz,NULL,NULL,iomodel_handle,"vz"); - IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements"); - IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements"); - - for (i=0;inumberofelements;i++){ - if(iomodel->my_elements[i]){ - - /*Create and add tria element to elements dataset: */ - elements->AddObject(new Penta(i+1,i,iomodel)); - - /*Create and add material property to materials dataset: */ - materials->AddObject(new Matice(i+1,i,iomodel,6)); - } - }//for (i=0;ielements); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->surface); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->drag_coefficient); - xfree((void**)&iomodel->drag_p); - xfree((void**)&iomodel->drag_q); - xfree((void**)&iomodel->elementoniceshelf); - xfree((void**)&iomodel->elementonbed); - xfree((void**)&iomodel->elementonsurface); - xfree((void**)&iomodel->elements_type); - xfree((void**)&iomodel->geothermalflux); - xfree((void**)&iomodel->rheology_n); - xfree((void**)&iomodel->rheology_B); - xfree((void**)&iomodel->elementonwater); - xfree((void**)&iomodel->pressure); - xfree((void**)&iomodel->temperature); - xfree((void**)&iomodel->vx); - xfree((void**)&iomodel->vy); - xfree((void**)&iomodel->vz); - xfree((void**)&iomodel->upperelements); - xfree((void**)&iomodel->lowerelements); - - /*Add new constrant material property tgo materials, at the end: */ - materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel));//put it at the end of the materials - - /*Create nodes and vertices: */ - if (iomodel->dim==3){ - IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); - IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); - } - IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); - IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); - IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); - IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); - IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); - IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); - IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); - IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); - IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); - - for (i=0;inumberofvertices;i++){ - - /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ - if(iomodel->my_vertices[i]){ - - /*Add vertex to vertices dataset: */ - vertices->AddObject(new Vertex(i+1,i,iomodel)); - - /*Add node to nodes dataset: */ - nodes->AddObject(new Node(i+1,i,iomodel)); - - } - } - - /*Clean fetched data: */ - xfree((void**)&iomodel->deadgrids); - xfree((void**)&iomodel->x); - xfree((void**)&iomodel->y); - xfree((void**)&iomodel->z); - xfree((void**)&iomodel->thickness); - xfree((void**)&iomodel->bed); - xfree((void**)&iomodel->gridonbed); - xfree((void**)&iomodel->gridonsurface); - xfree((void**)&iomodel->uppernodes); - xfree((void**)&iomodel->gridonicesheet); - xfree((void**)&iomodel->gridoniceshelf); - - /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these - * datasets, it will not be redone: */ - elements->Presort(); - nodes->Presort(); - vertices->Presort(); - materials->Presort(); - - cleanup_and_return: - - /*Assign output pointer: */ - *pelements=elements; - *pnodes=nodes; - *pvertices=vertices; - *pmaterials=materials; - -} Index: /issm/trunk/src/c/modules/ModelProcessorx/Thermal/CreateNodesThermal.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Thermal/CreateNodesThermal.cpp (revision 3998) +++ /issm/trunk/src/c/modules/ModelProcessorx/Thermal/CreateNodesThermal.cpp (revision 3998) @@ -0,0 +1,154 @@ +/* + * CreateElementsNodesAndMaterialsThermal.c: + */ + +#include "../../../DataSet/DataSet.h" +#include "../../../toolkits/toolkits.h" +#include "../../../io/io.h" +#include "../../../EnumDefinitions/EnumDefinitions.h" +#include "../../../objects/objects.h" +#include "../../../include/include.h" +#include "../../../shared/shared.h" +#include "../../MeshPartitionx/MeshPartitionx.h" +#include "../ModelProcessorx.h" + +void CreateElementsNodesAndMaterialsThermal(DataSet** pelements,DataSet** pnodes, DataSet** pvertices,DataSet** pmaterials, IoModel* iomodel,ConstDataHandle iomodel_handle){ + + /*Intermediary*/ + int i,j,k,n; + + /*DataSets: */ + DataSet* elements = NULL; + DataSet* nodes = NULL; + DataSet* vertices = NULL; + DataSet* materials = NULL; + + /*First create the elements, nodes and material properties: */ + elements = new DataSet(ElementsEnum); + nodes = new DataSet(NodesEnum); + vertices = new DataSet(VerticesEnum); + materials = new DataSet(MaterialsEnum); + + /*Partition elements and vertices and nodes: */ + Partitioning(&iomodel->my_elements, &iomodel->my_vertices, &iomodel->my_nodes, &iomodel->my_bordervertices, iomodel, iomodel_handle); + + /*Fetch data needed: */ + IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->drag_coefficient,NULL,NULL,iomodel_handle,"drag_coefficient"); + IoModelFetchData(&iomodel->drag_p,NULL,NULL,iomodel_handle,"drag_p"); + IoModelFetchData(&iomodel->drag_q,NULL,NULL,iomodel_handle,"drag_q"); + IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf"); + IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed"); + IoModelFetchData(&iomodel->elementonsurface,NULL,NULL,iomodel_handle,"elementonsurface"); + IoModelFetchData(&iomodel->elements_type,NULL,NULL,iomodel_handle,"elements_type"); + IoModelFetchData(&iomodel->geothermalflux,NULL,NULL,iomodel_handle,"geothermalflux"); + IoModelFetchData(&iomodel->rheology_B,NULL,NULL,iomodel_handle,"rheology_B"); + IoModelFetchData(&iomodel->rheology_n,NULL,NULL,iomodel_handle,"rheology_n"); + IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater"); + IoModelFetchData(&iomodel->pressure,NULL,NULL,iomodel_handle,"pressure"); + IoModelFetchData(&iomodel->temperature,NULL,NULL,iomodel_handle,"temperature"); + IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx"); + IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy"); + IoModelFetchData(&iomodel->vz,NULL,NULL,iomodel_handle,"vz"); + IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements"); + IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements"); + + for (i=0;inumberofelements;i++){ + if(iomodel->my_elements[i]){ + + /*Create and add tria element to elements dataset: */ + elements->AddObject(new Penta(i+1,i,iomodel)); + + /*Create and add material property to materials dataset: */ + materials->AddObject(new Matice(i+1,i,iomodel,6)); + } + }//for (i=0;ielements); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->surface); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->drag_coefficient); + xfree((void**)&iomodel->drag_p); + xfree((void**)&iomodel->drag_q); + xfree((void**)&iomodel->elementoniceshelf); + xfree((void**)&iomodel->elementonbed); + xfree((void**)&iomodel->elementonsurface); + xfree((void**)&iomodel->elements_type); + xfree((void**)&iomodel->geothermalflux); + xfree((void**)&iomodel->rheology_n); + xfree((void**)&iomodel->rheology_B); + xfree((void**)&iomodel->elementonwater); + xfree((void**)&iomodel->pressure); + xfree((void**)&iomodel->temperature); + xfree((void**)&iomodel->vx); + xfree((void**)&iomodel->vy); + xfree((void**)&iomodel->vz); + xfree((void**)&iomodel->upperelements); + xfree((void**)&iomodel->lowerelements); + + /*Add new constrant material property tgo materials, at the end: */ + materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel));//put it at the end of the materials + + /*Create nodes and vertices: */ + if (iomodel->dim==3){ + IoModelFetchData(&iomodel->deadgrids,NULL,NULL,iomodel_handle,"deadgrids"); + IoModelFetchData(&iomodel->uppernodes,NULL,NULL,iomodel_handle,"uppergrids"); + } + IoModelFetchData(&iomodel->x,NULL,NULL,iomodel_handle,"x"); + IoModelFetchData(&iomodel->y,NULL,NULL,iomodel_handle,"y"); + IoModelFetchData(&iomodel->z,NULL,NULL,iomodel_handle,"z"); + IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness"); + IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed"); + IoModelFetchData(&iomodel->gridonbed,NULL,NULL,iomodel_handle,"gridonbed"); + IoModelFetchData(&iomodel->gridonsurface,NULL,NULL,iomodel_handle,"gridonsurface"); + IoModelFetchData(&iomodel->gridonicesheet,NULL,NULL,iomodel_handle,"gridonicesheet"); + IoModelFetchData(&iomodel->gridoniceshelf,NULL,NULL,iomodel_handle,"gridoniceshelf"); + + for (i=0;inumberofvertices;i++){ + + /*vertices and nodes (same number, as we are running continuous galerkin formulation: */ + if(iomodel->my_vertices[i]){ + + /*Add vertex to vertices dataset: */ + vertices->AddObject(new Vertex(i+1,i,iomodel)); + + /*Add node to nodes dataset: */ + nodes->AddObject(new Node(i+1,i,iomodel)); + + } + } + + /*Clean fetched data: */ + xfree((void**)&iomodel->deadgrids); + xfree((void**)&iomodel->x); + xfree((void**)&iomodel->y); + xfree((void**)&iomodel->z); + xfree((void**)&iomodel->thickness); + xfree((void**)&iomodel->bed); + xfree((void**)&iomodel->gridonbed); + xfree((void**)&iomodel->gridonsurface); + xfree((void**)&iomodel->uppernodes); + xfree((void**)&iomodel->gridonicesheet); + xfree((void**)&iomodel->gridoniceshelf); + + /*All our datasets are already order by ids. Set presort flag so that later on, when sorting is requested on these + * datasets, it will not be redone: */ + elements->Presort(); + nodes->Presort(); + vertices->Presort(); + materials->Presort(); + + cleanup_and_return: + + /*Assign output pointer: */ + *pelements=elements; + *pnodes=nodes; + *pvertices=vertices; + *pmaterials=materials; + +}