Index: /issm/trunk/src/c/modules/ModelProcessorx/Balancevelocities/CreateConstraintsBalancevelocities.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Balancevelocities/CreateConstraintsBalancevelocities.cpp (revision 8822) +++ /issm/trunk/src/c/modules/ModelProcessorx/Balancevelocities/CreateConstraintsBalancevelocities.cpp (revision 8823) @@ -28,5 +28,6 @@ /*Fetch data: */ - IoModelFetchData(&iomodel->spcvelocity,NULL,NULL,iomodel_handle,"spcvelocity"); + IoModelFetchData(&iomodel->spcvx,NULL,NULL,iomodel_handle,"spcvx"); + IoModelFetchData(&iomodel->spcvy,NULL,NULL,iomodel_handle,"spcvy"); count=1; //matlab indexing @@ -36,9 +37,9 @@ if((iomodel->my_vertices[i])){ - if ((int)iomodel->spcvelocity[6*i+0] && (int)iomodel->spcvelocity[6*i+1]){ //spc if vx and vy are constrained + if (!isnan(iomodel->spcvx[i]) && !isnan(iomodel->spcvy[i])){ //spc if vx and vy are constrained /*This node needs to be spc'd: */ constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1, - pow( pow(*(iomodel->spcvelocity+6*i+4),2.0) + pow(*(iomodel->spcvelocity+6*i+5),2.0) ,0.5),BalancevelocitiesAnalysisEnum)); + pow( pow(iomodel->spcvx[i],2.0) + pow(iomodel->spcvy[i],2.0) ,0.5),BalancevelocitiesAnalysisEnum)); count++; } @@ -47,5 +48,6 @@ /*Free data: */ - xfree((void**)&iomodel->spcvelocity); + xfree((void**)&iomodel->spcvx); + xfree((void**)&iomodel->spcvy); /*Assign output pointer: */ Index: /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHoriz/CreateConstraintsDiagnosticHoriz.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHoriz/CreateConstraintsDiagnosticHoriz.cpp (revision 8822) +++ /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHoriz/CreateConstraintsDiagnosticHoriz.cpp (revision 8823) @@ -32,5 +32,7 @@ /*Spcs: fetch data: */ - IoModelFetchData(&iomodel->spcvelocity,NULL,NULL,iomodel_handle,"spcvelocity"); + IoModelFetchData(&iomodel->spcvx,NULL,NULL,iomodel_handle,"spcvx"); + IoModelFetchData(&iomodel->spcvy,NULL,NULL,iomodel_handle,"spcvy"); + IoModelFetchData(&iomodel->spcvz,NULL,NULL,iomodel_handle,"spcvz"); IoModelFetchData(&iomodel->nodeonhutter,NULL,NULL,iomodel_handle,"nodeonhutter"); IoModelFetchData(&iomodel->nodeonmacayeal,NULL,NULL,iomodel_handle,"nodeonmacayeal"); @@ -56,10 +58,10 @@ constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,0,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. count++; - if ((int)iomodel->spcvelocity[6*i+0]){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,*(iomodel->spcvelocity+6*i+3)/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. - count++; - } - if ((int)iomodel->spcvelocity[6*i+1]){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,*(iomodel->spcvelocity+6*i+4)/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. + if (!isnan(iomodel->spcvx[i])){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,iomodel->spcvx[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. + count++; + } + if (!isnan(iomodel->spcvy[i])){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,iomodel->spcvy[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. count++; } @@ -71,10 +73,10 @@ constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,0,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. count++; - if ((int)iomodel->spcvelocity[6*i+0]){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,*(iomodel->spcvelocity+6*i+3)/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. - count++; - } - if ((int)iomodel->spcvelocity[6*i+1]){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,*(iomodel->spcvelocity+6*i+4)/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. + if (!isnan(iomodel->spcvx[i])){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,iomodel->spcvx[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. + count++; + } + if (!isnan(iomodel->spcvy[i])){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,iomodel->spcvy[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. count++; } @@ -93,10 +95,10 @@ constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,5,0,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. count++; - if ((int)iomodel->spcvelocity[6*i+0]){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,*(iomodel->spcvelocity+6*i+3)/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. - count++; - } - if ((int)iomodel->spcvelocity[6*i+1]){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,*(iomodel->spcvelocity+6*i+4)/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. + if (!isnan(iomodel->spcvx[i])){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,iomodel->spcvx[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. + count++; + } + if (!isnan(iomodel->spcvy[i])){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,iomodel->spcvy[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. count++; } @@ -108,14 +110,14 @@ constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,0,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. count++; - if ((int)iomodel->spcvelocity[6*i+0]){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,*(iomodel->spcvelocity+6*i+3)/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. - count++; - } - if ((int)iomodel->spcvelocity[6*i+1]){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,*(iomodel->spcvelocity+6*i+4)/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. - count++; - } - if ((int)iomodel->spcvelocity[6*i+2]){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,5,*(iomodel->spcvelocity+6*i+5)/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. + if (!isnan(iomodel->spcvx[i])){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,iomodel->spcvx[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. + count++; + } + if (!isnan(iomodel->spcvy[i])){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,iomodel->spcvy[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. + count++; + } + if (!isnan(iomodel->spcvz[i])){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,5,iomodel->spcvz[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. count++; } @@ -133,10 +135,10 @@ constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,5,0,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. count++; - if ((int)iomodel->spcvelocity[6*i+0]){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,*(iomodel->spcvelocity+6*i+3)/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. - count++; - } - if ((int)iomodel->spcvelocity[6*i+1]){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,*(iomodel->spcvelocity+6*i+4)/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. + if (!isnan(iomodel->spcvx[i])){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,iomodel->spcvx[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. + count++; + } + if (!isnan(iomodel->spcvy[i])){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,iomodel->spcvy[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. count++; } @@ -148,14 +150,14 @@ constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,0,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. count++; - if ((int)iomodel->spcvelocity[6*i+0]){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,*(iomodel->spcvelocity+6*i+3)/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. - count++; - } - if ((int)iomodel->spcvelocity[6*i+1]){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,*(iomodel->spcvelocity+6*i+4)/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. - count++; - } - if ((int)iomodel->spcvelocity[6*i+2]){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,5,*(iomodel->spcvelocity+6*i+5)/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. + if (!isnan(iomodel->spcvx[i])){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,iomodel->spcvx[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. + count++; + } + if (!isnan(iomodel->spcvy[i])){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,iomodel->spcvy[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. + count++; + } + if (!isnan(iomodel->spcvz[i])){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,5,iomodel->spcvz[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. count++; } @@ -165,15 +167,15 @@ /*Now add the regular spcs*/ else{ - if ((int)iomodel->spcvelocity[6*i+0] || (int)iomodel->nodeonhutter[i]){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,*(iomodel->spcvelocity+6*i+3)/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. + if (!isnan(iomodel->spcvx[i]) || (int)iomodel->nodeonhutter[i]){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,iomodel->spcvx[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. count++; } - if ((int)iomodel->spcvelocity[6*i+1] || (int)iomodel->nodeonhutter[i]){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,*(iomodel->spcvelocity+6*i+4)/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 2 to vy - count++; - } - if ((int)iomodel->spcvelocity[6*i+2] && ((int)iomodel->vertices_type[i]==StokesApproximationEnum || ((int)iomodel->vertices_type[i]==NoneApproximationEnum))){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,*(iomodel->spcvelocity+6*i+5)/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 2 to vy + if (!isnan(iomodel->spcvy[i]) || (int)iomodel->nodeonhutter[i]){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,iomodel->spcvy[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 2 to vy + count++; + } + if (!isnan(iomodel->spcvz[i]) && ((int)iomodel->vertices_type[i]==StokesApproximationEnum || ((int)iomodel->vertices_type[i]==NoneApproximationEnum))){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,iomodel->spcvz[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 2 to vy count++; } @@ -187,5 +189,7 @@ /*Free data: */ - xfree((void**)&iomodel->spcvelocity); + xfree((void**)&iomodel->spcvx); + xfree((void**)&iomodel->spcvy); + xfree((void**)&iomodel->spcvz); xfree((void**)&iomodel->nodeonhutter); xfree((void**)&iomodel->nodeonmacayeal); Index: /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHoriz/CreateLoadsDiagnosticHoriz.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHoriz/CreateLoadsDiagnosticHoriz.cpp (revision 8822) +++ /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHoriz/CreateLoadsDiagnosticHoriz.cpp (revision 8823) @@ -107,5 +107,7 @@ IoModelFetchData(&iomodel->nodeonstokes,NULL,NULL,iomodel_handle,"nodeonstokes"); IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements"); - IoModelFetchData(&iomodel->spcvelocity,NULL,NULL,iomodel_handle,"spcvelocity"); + IoModelFetchData(&iomodel->spcvx,NULL,NULL,iomodel_handle,"spcvx"); + IoModelFetchData(&iomodel->spcvy,NULL,NULL,iomodel_handle,"spcvy"); + IoModelFetchData(&iomodel->spcvz,NULL,NULL,iomodel_handle,"spcvz"); IoModelFetchData(&iomodel->vertices_type,NULL,NULL,iomodel_handle,"vertices_type"); CreateSingleNodeToElementConnectivity(iomodel); @@ -128,5 +130,7 @@ xfree((void**)&iomodel->nodeonicesheet); xfree((void**)&iomodel->elements); - xfree((void**)&iomodel->spcvelocity); + xfree((void**)&iomodel->spcvx); + xfree((void**)&iomodel->spcvy); + xfree((void**)&iomodel->spcvz); xfree((void**)&iomodel->vertices_type); xfree((void**)&iomodel->singlenodetoelementconnectivity); Index: /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHutter/CreateConstraintsDiagnosticHutter.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHutter/CreateConstraintsDiagnosticHutter.cpp (revision 8822) +++ /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHutter/CreateConstraintsDiagnosticHutter.cpp (revision 8823) @@ -31,5 +31,6 @@ /*Fetch data: */ - IoModelFetchData(&iomodel->spcvelocity,NULL,NULL,iomodel_handle,"spcvelocity"); + IoModelFetchData(&iomodel->spcvx,NULL,NULL,iomodel_handle,"spcvx"); + IoModelFetchData(&iomodel->spcvy,NULL,NULL,iomodel_handle,"spcvy"); IoModelFetchData(&iomodel->nodeonhutter,NULL,NULL,iomodel_handle,"nodeonhutter"); @@ -50,11 +51,11 @@ } else{ - if ((int)iomodel->spcvelocity[6*i+0]){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,*(iomodel->spcvelocity+6*i+3)/iomodel->yts,DiagnosticHutterAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. + if (!isnan(iomodel->spcvx[i])){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,iomodel->spcvx[i]/iomodel->yts,DiagnosticHutterAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. count++; } - if ((int)iomodel->spcvelocity[6*i+1]){ - constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,*(iomodel->spcvelocity+6*i+4)/iomodel->yts,DiagnosticHutterAnalysisEnum)); //add count'th spc, on node i+1, setting dof 2 to vy + if (!isnan(iomodel->spcvy[i])){ + constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,iomodel->spcvy[i]/iomodel->yts,DiagnosticHutterAnalysisEnum)); //add count'th spc, on node i+1, setting dof 2 to vy count++; } @@ -65,5 +66,6 @@ /*Free data: */ xfree((void**)&iomodel->nodeonhutter); - xfree((void**)&iomodel->spcvelocity); + xfree((void**)&iomodel->spcvx); + xfree((void**)&iomodel->spcvy); cleanup_and_return: Index: /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticVert/CreateConstraintsDiagnosticVert.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticVert/CreateConstraintsDiagnosticVert.cpp (revision 8822) +++ /issm/trunk/src/c/modules/ModelProcessorx/DiagnosticVert/CreateConstraintsDiagnosticVert.cpp (revision 8823) @@ -30,5 +30,5 @@ /*Fetch data: */ - IoModelFetchData(&iomodel->spcvelocity,NULL,NULL,iomodel_handle,"spcvelocity"); + IoModelFetchData(&iomodel->spcvz,NULL,NULL,iomodel_handle,"spcvz"); IoModelFetchData(&iomodel->nodeonstokes,NULL,NULL,iomodel_handle,"nodeonstokes"); @@ -46,7 +46,7 @@ count++; } - else if ((int)iomodel->spcvelocity[6*i+2]){ + else if (!isnan(iomodel->spcvz[i])){ constraints->AddObject(new Spc(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1, - *(iomodel->spcvelocity+6*i+5)/iomodel->yts,DiagnosticVertAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. + iomodel->spcvz[i]/iomodel->yts,DiagnosticVertAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx. count++; @@ -56,5 +56,5 @@ /*Free data: */ - xfree((void**)&iomodel->spcvelocity); + xfree((void**)&iomodel->spcvz); xfree((void**)&iomodel->nodeonstokes); Index: /issm/trunk/src/c/objects/IoModel.cpp =================================================================== --- /issm/trunk/src/c/objects/IoModel.cpp (revision 8822) +++ /issm/trunk/src/c/objects/IoModel.cpp (revision 8823) @@ -79,5 +79,7 @@ xfree((void**)&this->nodeonwater); xfree((void**)&this->pressureload); - xfree((void**)&this->spcvelocity); + xfree((void**)&this->spcvx); + xfree((void**)&this->spcvy); + xfree((void**)&this->spcvz); xfree((void**)&this->spcthickness); xfree((void**)&this->spcwatercolumn); @@ -325,5 +327,7 @@ this->numberofpressureloads=0; this->pressureload=NULL; - this-> spcvelocity=NULL; + this-> spcvx=NULL; + this-> spcvy=NULL; + this-> spcvz=NULL; this-> spctemperature=NULL; this-> spcthickness=NULL; Index: /issm/trunk/src/c/objects/IoModel.h =================================================================== --- /issm/trunk/src/c/objects/IoModel.h (revision 8822) +++ /issm/trunk/src/c/objects/IoModel.h (revision 8823) @@ -101,5 +101,7 @@ int numberofpressureloads; double* pressureload; - double* spcvelocity; + double* spcvx; + double* spcvy; + double* spcvz; double* spctemperature; double* spcthickness; Index: /issm/trunk/src/m/classes/version/7.6/model.m =================================================================== --- /issm/trunk/src/m/classes/version/7.6/model.m (revision 8822) +++ /issm/trunk/src/m/classes/version/7.6/model.m (revision 8823) @@ -161,5 +161,7 @@ nodeonboundary=NaN; pressureload=NaN; - spcvelocity=NaN; + spcvx=NaN; + spcvy=NaN; + spcvz=NaN; spctemperature=NaN; spcthickness=NaN; Index: /issm/trunk/src/m/classes/version/7.7/model.m =================================================================== --- /issm/trunk/src/m/classes/version/7.7/model.m (revision 8822) +++ /issm/trunk/src/m/classes/version/7.7/model.m (revision 8823) @@ -161,5 +161,7 @@ nodeonboundary=NaN; pressureload=NaN; - spcvelocity=NaN; + spcvx=NaN; + spcvy=NaN; + spcvz=NaN; spctemperature=NaN; spcthickness=NaN; Index: /issm/trunk/src/m/model/BasinConstrain.m =================================================================== --- /issm/trunk/src/m/model/BasinConstrain.m (revision 8822) +++ /issm/trunk/src/m/model/BasinConstrain.m (revision 8823) @@ -48,7 +48,6 @@ %all elements outside the constraint domain are equivalent to water. all nodes outside are spc'd. -md.spcvelocity(nodenotondomain,1:2)=1; -md.spcvelocity(nodenotondomain,4)=md.vx_obs(nodenotondomain); -md.spcvelocity(nodenotondomain,5)=md.vy_obs(nodenotondomain); +md.spcvx(nodenotondomain)=md.vx_obs(nodenotondomain); +md.spcvy(nodenotondomain)=md.vy_obs(nodenotondomain); md.elementonwater(elementnotondomain)=1; @@ -57,9 +56,8 @@ numpos=unique(md.elements(pos,:)); nodes=setdiff(1:1:md.numberofnodes,numpos); -md.spcvelocity(nodes,1:2)=1; -md.spcvelocity(nodes,4)=md.vx_obs(nodes); -md.spcvelocity(nodes,5)=md.vy_obs(nodes); +md.spcvx(nodes)=md.vx_obs(nodes); +md.spcvy(nodes)=md.vy_obs(nodes); %make sure icefronts that are completely spc'd are taken out: -free_segments=find(sum(md.spcvelocity(md.pressureload(:,1:2),1:2),2)~=2); +free_segments=find((~isnan(md.spcvx(md.pressureload(:,1:2))) + ~isnan(md.spcvy(md.pressureload(:,1:2))))~=2); md.pressureload=md.pressureload(free_segments,:); Index: /issm/trunk/src/m/model/BasinConstrain2.m =================================================================== --- /issm/trunk/src/m/model/BasinConstrain2.m (revision 8822) +++ /issm/trunk/src/m/model/BasinConstrain2.m (revision 8823) @@ -48,7 +48,6 @@ %all elements outside the constraint domain are equivalent to water. all nodes outside are spc'd. -md.spcvelocity(nodenotondomain,1:2)=1; -md.spcvelocity(nodenotondomain,4)=md.vx_obs(nodenotondomain); -md.spcvelocity(nodenotondomain,5)=md.vy_obs(nodenotondomain); +md.spcvx(nodenotondomain)=md.vx_obs(nodenotondomain); +md.spcvy(nodenotondomain)=md.vy_obs(nodenotondomain); md.elementonwater(elementnotondomain)=1; @@ -57,10 +56,9 @@ numpos=unique(md.elements(pos,:)); nodes=setdiff(1:1:md.numberofnodes,numpos); -md.spcvelocity(nodes,1:2)=1; -md.spcvelocity(nodes,4)=md.vx_obs(nodes); -md.spcvelocity(nodes,5)=md.vy_obs(nodes); +md.spcvx(nodes)=md.vx_obs(nodes); +md.spcvy(nodes)=md.vy_obs(nodes); %make sure icefronts that are completely spc'd are taken out: -free_segments=find(sum(md.spcvelocity(md.pressureload(:,1:2),1:2),2)~=2); +free_segments=find((~isnan(md.spcvx(md.pressureload(:,1:2))) + ~isnan(md.spcvy(md.pressureload(:,1:2))))~=2); md.pressureload=md.pressureload(free_segments,:); Index: /issm/trunk/src/m/model/BasinConstrainShelf.m =================================================================== --- /issm/trunk/src/m/model/BasinConstrainShelf.m (revision 8822) +++ /issm/trunk/src/m/model/BasinConstrainShelf.m (revision 8823) @@ -48,7 +48,6 @@ %all elements outside the constraint domain are equivalent to water. all nodes outside are spc'd. -md.spcvelocity(nodenotondomain,1:2)=1; -md.spcvelocity(nodenotondomain,4)=md.vx_obs(nodenotondomain); -md.spcvelocity(nodenotondomain,5)=md.vy_obs(nodenotondomain); +md.spcvx(nodenotondomain)=md.vx_obs(nodenotondomain); +md.spcvy(nodenotondomain)=md.vy_obs(nodenotondomain); md.elementonwater(elementnotondomain)=1; @@ -57,23 +56,20 @@ numpos=unique(md.elements(pos,:)); nodes=setdiff(1:1:md.numberofnodes,numpos); -md.spcvelocity(nodes,1:2)=1; -md.spcvelocity(nodes,4)=md.vx_obs(nodes); -md.spcvelocity(nodes,5)=md.vy_obs(nodes); +md.spcvx(nodes)=md.vx_obs(nodes); +md.spcvy(nodes)=md.vy_obs(nodes); %make sure any node with NaN velocity is spc'd: +%we spc to the smoothed value, so that control methods don't go berserk trying to figure out what reaction force to apply for the spc to stand. pos=find(isnan(md.vel_obs_raw)); -md.spcvelocity(pos,1:2)=1; -%we spc to the smoothed value, so that control methods don't go berserk trying to figure out what reaction force to apply for the spc to stand. -md.spcvelocity(pos,4)=md.vx_obs(pos); -md.spcvelocity(pos,5)=md.vy_obs(pos); +md.spcvx(pos)=md.vx_obs(pos); +md.spcvy(pos)=md.vy_obs(pos); %iceshelves: any node on icesheet is spc'd pos=find(md.nodeonicesheet); -md.spcvelocity(pos,1:2)=1; -md.spcvelocity(pos,4)=md.vx_obs(pos); -md.spcvelocity(pos,5)=md.vy_obs(pos); +md.spcvx(pos)=md.vx_obs(pos); +md.spcvy(pos)=md.vy_obs(pos); %make sure icefronts that are completely spc'd are taken out: -free_segments=find(sum(md.spcvelocity(md.pressureload(:,1:2),1:2),2)~=2); +free_segments=find((~isnan(md.spcvx(md.pressureload(:,1:2))) + ~isnan(md.spcvy(md.pressureload(:,1:2))) )~=2); md.pressureload=md.pressureload(free_segments,:); Index: /issm/trunk/src/m/model/collapse.m =================================================================== --- /issm/trunk/src/m/model/collapse.m (revision 8822) +++ /issm/trunk/src/m/model/collapse.m (revision 8823) @@ -61,5 +61,7 @@ %boundary conditions -md.spcvelocity=project2d(md,md.spcvelocity,md.numlayers); +md.spcvx=project2d(md,md.spcvx,md.numlayers); +md.spcvy=project2d(md,md.spcvy,md.numlayers); +md.spcvz=project2d(md,md.spcvz,md.numlayers); md.spcthickness=project2d(md,md.spcthickness,md.numlayers); md.spctemperature=project2d(md,md.spctemperature,md.numlayers); Index: /issm/trunk/src/m/model/display/displaybc.m =================================================================== --- /issm/trunk/src/m/model/display/displaybc.m (revision 8822) +++ /issm/trunk/src/m/model/display/displaybc.m (revision 8823) @@ -19,5 +19,7 @@ disp(sprintf('\n diagnostic:')); -fielddisplay(md,'spcvelocity','constraints flag list (first 3 columns) and values [m/yr] (last 3 columns)'); +fielddisplay(md,'spcvx','x-axis velocity constraint (NaN means no constraint)'); +fielddisplay(md,'spcvy','y-axis velocity constraint (NaN means no constraint)'); +fielddisplay(md,'spcvz','z-axis velocity constraint (NaN means no constraint)'); fielddisplay(md,'pressureload','segments on ice front list'); Index: /issm/trunk/src/m/model/display/displaydiagnostic.m =================================================================== --- /issm/trunk/src/m/model/display/displaydiagnostic.m (revision 8822) +++ /issm/trunk/src/m/model/display/displaydiagnostic.m (revision 8823) @@ -19,5 +19,7 @@ disp(sprintf('\n boundary conditions:')); -fielddisplay(md,'spcvelocity','constraints flag list (first 3 columns) and values [m/yr] (last 3 columns)'); +fielddisplay(md,'spcvx','x-axis velocity constraint (NaN means no constraint)'); +fielddisplay(md,'spcvy','y-axis velocity constraint (NaN means no constraint)'); +fielddisplay(md,'spcvz','z-axis velocity constraint (NaN means no constraint)'); fielddisplay(md,'pressureload','segments on ice front list'); Index: /issm/trunk/src/m/model/extrude.m =================================================================== --- /issm/trunk/src/m/model/extrude.m (revision 8822) +++ /issm/trunk/src/m/model/extrude.m (revision 8823) @@ -201,5 +201,7 @@ %boundary conditions -md.spcvelocity=project3d(md,md.spcvelocity,'node'); +md.spcvx=project3d(md,md.spcvx,'node'); +md.spcvy=project3d(md,md.spcvy,'node'); +md.spcvz=project3d(md,md.spcvz,'node'); md.spctemperature=project3d(md,md.spctemperature,'node',md.numlayers); md.spcthickness=project3d(md,md.spcthickness,'node'); Index: /issm/trunk/src/m/model/ismodelselfconsistent.m =================================================================== --- /issm/trunk/src/m/model/ismodelselfconsistent.m (revision 8822) +++ /issm/trunk/src/m/model/ismodelselfconsistent.m (revision 8823) @@ -99,5 +99,5 @@ %}}} %OTHER SIZES {{{1 -fields={'spcvelocity','diagnostic_ref'}; +fields={'diagnostic_ref'}; checksize(md,fields,[md.numberofnodes 6]); %}}} @@ -385,5 +385,5 @@ % {{{2 %SINGULAR - if ~any(sum(md.spcvelocity(:,1:2),2)==2), + if ~any((~isnan(md.spcvx)+~isnan(md.spcvy))==2), message(['model not consistent: model ' md.name ' is not well posed (singular). You need at least one node with fixed velocity!']) end @@ -603,6 +603,10 @@ %SPC - if any(md.spcvelocity(find(md.nodeonboundary),[1:2])~=1), - message(['model not consistent: model ' md.name ' should have all the nodes on boundary constrained in field spcvelocity']); + if any((~isnan(md.spcvx(find(md.nodeonboundary))))~=1), + message(['model not consistent: model ' md.name ' should have all the nodes on boundary constrained in field spcvx']); + end + + if any((~isnan(md.spcvy(find(md.nodeonboundary))))~=1), + message(['model not consistent: model ' md.name ' should have all the nodes on boundary constrained in field spcvy']); end %}}} Index: /issm/trunk/src/m/model/marshall.m =================================================================== --- /issm/trunk/src/m/model/marshall.m (revision 8822) +++ /issm/trunk/src/m/model/marshall.m (revision 8823) @@ -77,5 +77,7 @@ WriteData(fid,md.nodeonwater,'Mat','nodeonwater'); -WriteData(fid,md.spcvelocity,'Mat','spcvelocity'); +WriteData(fid,md.spcvx,'Mat','spcvx'); +WriteData(fid,md.spcvy,'Mat','spcvy'); +WriteData(fid,md.spcvz,'Mat','spcvz'); WriteData(fid,md.spctemperature,'Mat','spctemperature'); WriteData(fid,md.spcthickness,'Mat','spcthickness'); Index: /issm/trunk/src/m/model/modelextract.m =================================================================== --- /issm/trunk/src/m/model/modelextract.m (revision 8822) +++ /issm/trunk/src/m/model/modelextract.m (revision 8823) @@ -195,11 +195,11 @@ nodestoflag1=intersect(orphans_node,pos_node); nodestoflag2=Pnode(nodestoflag1); - if ~isnan(md1.spcvelocity), - md2.spcvelocity(nodestoflag2,1:3)=1; + if ~isnan(md1.spcvx) && ~isnan(md1.spcvy) && ~isnan(md1.spcvz), if ~isnan(md1.vx_obs) & ~isnan(md1.vy_obs) - md2.spcvelocity(nodestoflag2,4)=md2.vx_obs(nodestoflag2); - md2.spcvelocity(nodestoflag2,5)=md2.vy_obs(nodestoflag2); + md2.spcvx(nodestoflag2)=md2.vx_obs(nodestoflag2); + md2.spcvy(nodestoflag2)=md2.vy_obs(nodestoflag2); else - md2.spcvelocity(nodestoflag2,4:5)=zeros(length(nodestoflag2),2); + md2.spcvx(nodestoflag2)=NaN; + md2.spcvy(nodestoflag2)=NaN; disp(' ') disp('!! modelextract warning: spc values should be checked !!') Index: /issm/trunk/src/m/model/plot/plot_BC.m =================================================================== --- /issm/trunk/src/m/model/plot/plot_BC.m (revision 8822) +++ /issm/trunk/src/m/model/plot/plot_BC.m (revision 8823) @@ -7,7 +7,7 @@ %plot dirichlets -h1=plot(md.x(find(md.spcvelocity(:,1))),md.y(find(md.spcvelocity(:,1))),'ro','MarkerSize',14,'MarkerFaceColor','r'); -h2=plot(md.x(find(md.spcvelocity(:,2))),md.y(find(md.spcvelocity(:,2))),'bo','MarkerSize',10,'MarkerFaceColor','b'); -h3=plot(md.x(find(md.spcvelocity(:,3))),md.y(find(md.spcvelocity(:,3))),'yo','MarkerSize',6 ,'MarkerFaceColor','y'); +h1=plot(md.x(find(~isnan(md.spcvx))),md.y(find(~isnan(md.spcvx))),'ro','MarkerSize',14,'MarkerFaceColor','r'); +h2=plot(md.x(find(~isnan(md.spcvy))),md.y(find(~isnan(md.spcvy))),'bo','MarkerSize',10,'MarkerFaceColor','b'); +h3=plot(md.x(find(~isnan(md.spcvz))),md.y(find(~isnan(md.spcvz))),'yo','MarkerSize',6 ,'MarkerFaceColor','y'); %update legend Index: /issm/trunk/src/m/model/setelementstype.m =================================================================== --- /issm/trunk/src/m/model/setelementstype.m (revision 8822) +++ /issm/trunk/src/m/model/setelementstype.m (revision 8823) @@ -97,5 +97,5 @@ nodeonstokes=zeros(md.numberofnodes,1); nodeonstokes(md.elements(find(stokesflag),:))=1; - fullspcnodes=double(sum(md.spcvelocity(:,1:3),2)==3 | (nodeonpattyn & nodeonstokes)); %find all the nodes on the boundary of the domain without icefront + fullspcnodes=double((~isnan(md.spcvx)+~isnan(md.spcvy)+~isnan(md.spcvz))==3 | (nodeonpattyn & nodeonstokes)); %find all the nodes on the boundary of the domain without icefront fullspcelems=double(sum(fullspcnodes(md.elements),2)==6); %find all the nodes on the boundary of the domain without icefront stokesflag(find(fullspcelems))=0; Index: /issm/trunk/src/m/utils/Analysis/resetspcs.m =================================================================== --- /issm/trunk/src/m/utils/Analysis/resetspcs.m (revision 8822) +++ /issm/trunk/src/m/utils/Analysis/resetspcs.m (revision 8823) @@ -3,6 +3,10 @@ resolution=20000; x_ellsworth=-1.24*10^6; y_ellsworth=1*10^5; pos=find((md.x>(x_ellsworth-resolution)) & (md.x<(x_ellsworth+resolution)) & (md.y>(y_ellsworth-resolution)) & (md.y<(y_ellsworth+resolution))); - md.spcvelocity=zeros(md.numberofnodes,6); - md.spcvelocity(pos,1:3)=1; + md.spcvx=NaN*ones(md.numberofnodes); + md.spcvy=NaN*ones(md.numberofnodes); + md.spcvz=NaN*ones(md.numberofnodes); + md.spcvx(pos)=0; + md.spcvy(pos)=0; + md.spcvz(pos)=0; %recompute with spcs reinitialized Index: /issm/trunk/src/m/utils/BC/SetIceSheetBC.m =================================================================== --- /issm/trunk/src/m/utils/BC/SetIceSheetBC.m (revision 8822) +++ /issm/trunk/src/m/utils/BC/SetIceSheetBC.m (revision 8823) @@ -9,6 +9,10 @@ %node on Dirichlet pos=find(md.nodeonboundary); -md.spcvelocity=zeros(md.numberofnodes,6); -md.spcvelocity(pos,1:3)=1; +md.spcvx=NaN*ones(md.numberofnodes,1); +md.spcvy=NaN*ones(md.numberofnodes,1); +md.spcvz=NaN*ones(md.numberofnodes,1); +md.spcvx(pos)=0; +md.spcvy(pos)=0; +md.spcvz(pos)=0; md.diagnostic_ref=NaN*ones(md.numberofnodes,6); @@ -16,5 +20,6 @@ if (length(md.vx_obs)==md.numberofnodes & length(md.vy_obs)==md.numberofnodes) disp(' boundary conditions for diagnostic model: spc set as observed velocities'); - md.spcvelocity(:,4:5)=[md.vx_obs md.vy_obs]; %vz is zero + md.spcvx=md.vx_obs(pos); + md.spcvy=md.vy_obs(pos); else disp(' boundary conditions for diagnostic model: spc set as zero'); Index: /issm/trunk/src/m/utils/BC/SetIceShelfBC.m =================================================================== --- /issm/trunk/src/m/utils/BC/SetIceShelfBC.m (revision 8822) +++ /issm/trunk/src/m/utils/BC/SetIceShelfBC.m (revision 8823) @@ -28,6 +28,10 @@ end pos=find(md.nodeonboundary & ~nodeonicefront); -md.spcvelocity=zeros(md.numberofnodes,6); -md.spcvelocity(pos,1:3)=1; +md.spcvx=NaN*ones(md.numberofnodes,1); +md.spcvy=NaN*ones(md.numberofnodes,1); +md.spcvz=NaN*ones(md.numberofnodes,1); +md.spcvx(pos)=0; +md.spcvy(pos)=0; +md.spcvz(pos)=0; md.diagnostic_ref=NaN*ones(md.numberofnodes,6); @@ -35,5 +39,6 @@ if (length(md.vx_obs)==md.numberofnodes & length(md.vy_obs)==md.numberofnodes) disp(' boundary conditions for diagnostic model: spc set as observed velocities'); - md.spcvelocity(pos,4:5)=[md.vx_obs(pos) md.vy_obs(pos)]; %zeros for vz + md.spcvx(pos)=md.vx_obs(pos); + md.spcvy(pos)=md.vy_obs(pos); else disp(' boundary conditions for diagnostic model: spc set as zero'); Index: /issm/trunk/src/m/utils/BC/SetMarineIceSheetBC.m =================================================================== --- /issm/trunk/src/m/utils/BC/SetMarineIceSheetBC.m (revision 8822) +++ /issm/trunk/src/m/utils/BC/SetMarineIceSheetBC.m (revision 8823) @@ -35,6 +35,10 @@ warning('SetMarineIceSheetBC warning: ice front all around the glacier, no dirichlet found. Dirichlet must be added manually') end -md.spcvelocity=zeros(md.numberofnodes,6); -md.spcvelocity(pos,1:3)=1; +md.spcvx=NaN*ones(md.numberofnodes,1); +md.spcvy=NaN*ones(md.numberofnodes,1); +md.spcvz=NaN*ones(md.numberofnodes,1); +md.spcvx(pos)=0; +md.spcvy(pos)=0; +md.spcvz(pos)=0; md.diagnostic_ref=NaN*ones(md.numberofnodes,6); @@ -42,5 +46,6 @@ if (length(md.vx_obs)==md.numberofnodes & length(md.vy_obs)==md.numberofnodes) disp(' boundary conditions for diagnostic model: spc set as observed velocities'); - md.spcvelocity(pos,4:5)=[md.vx_obs(pos) md.vy_obs(pos)]; %zeros for vz + md.spcvx(pos)=md.vx_obs(pos); + md.spcvy(pos)=md.vy_obs(pos); else disp(' boundary conditions for diagnostic model: spc set as zero');