Index: /issm/trunk-jpl/src/c/Makefile.am =================================================================== --- /issm/trunk-jpl/src/c/Makefile.am (revision 13874) +++ /issm/trunk-jpl/src/c/Makefile.am (revision 13875) @@ -481,12 +481,4 @@ ./solutions/balancethickness_core.cpp #}}} -#Responses sources {{{ -responses_sources = ./modules/TotalSmbx/TotalSmbx.h\ - ./modules/TotalSmbx/TotalSmbx.cpp\ - ./modules/IceVolumex/IceVolumex.h\ - ./modules/IceVolumex/IceVolumex.cpp\ - ./modules/ElementResponsex/ElementResponsex.h\ - ./modules/ElementResponsex/ElementResponsex.cpp -#}}} #Slope sources {{{ slope_sources = ./modules/ModelProcessorx/BedSlope/UpdateElementsBedSlope.cpp\ @@ -863,8 +855,4 @@ endif -if RESPONSES -issm_sources += $(responses_sources) -endif - if SLOPE issm_sources += $(slope_sources) Index: /issm/trunk-jpl/src/c/classes/FemModel.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/FemModel.cpp (revision 13874) +++ /issm/trunk-jpl/src/c/classes/FemModel.cpp (revision 13875) @@ -404,5 +404,5 @@ #ifdef _HAVE_RESPONSES_ - case IceVolumeEnum: IceVolumex( responses, elements,nodes, vertices, loads, materials, parameters,process_units); break; + case IceVolumeEnum: this->IceVolumex(responses,process_units); break; case MinVelEnum: this->MinVelx(responses,process_units); break; case MaxVelEnum: this->MaxVelx( responses,process_units); break; @@ -416,19 +416,19 @@ case MaxVzEnum: this->MaxVzx( responses,process_units); break; case MaxAbsVzEnum: this->MaxAbsVzx( responses,process_units); break; - case MassFluxEnum: this->MassFluxx( responses,process_units); break; - case SurfaceAbsVelMisfitEnum: SurfaceAbsVelMisfitx( responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break; - case SurfaceRelVelMisfitEnum: SurfaceRelVelMisfitx( responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break; - case SurfaceLogVelMisfitEnum: SurfaceLogVelMisfitx( responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break; - case SurfaceLogVxVyMisfitEnum: SurfaceLogVxVyMisfitx( responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break; + case MassFluxEnum: this->MassFluxx( responses,process_units); break; + case SurfaceAbsVelMisfitEnum: SurfaceAbsVelMisfitx(responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break; + case SurfaceRelVelMisfitEnum: SurfaceRelVelMisfitx(responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break; + case SurfaceLogVelMisfitEnum: SurfaceLogVelMisfitx(responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break; + case SurfaceLogVxVyMisfitEnum: SurfaceLogVxVyMisfitx(responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break; case SurfaceAverageVelMisfitEnum:SurfaceAverageVelMisfitx( responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break; - case ThicknessAbsMisfitEnum: ThicknessAbsMisfitx( responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break; - case ThicknessAbsGradientEnum: this->ThicknessAbsGradientx( responses, process_units,weight_index); break; - case ThicknessAlongGradientEnum: ThicknessAlongGradientx( responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break; - case ThicknessAcrossGradientEnum:ThicknessAcrossGradientx( responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break; - case TotalSmbEnum: TotalSmbx( responses, elements,nodes, vertices, loads, materials, parameters,process_units); break; - case RheologyBbarAbsGradientEnum:RheologyBbarAbsGradientx( responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break; + case ThicknessAbsMisfitEnum: ThicknessAbsMisfitx(responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break; + case ThicknessAbsGradientEnum: this->ThicknessAbsGradientx(responses, process_units,weight_index); break; + case ThicknessAlongGradientEnum: ThicknessAlongGradientx(responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break; + case ThicknessAcrossGradientEnum:ThicknessAcrossGradientx(responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break; + case TotalSmbEnum: this->TotalSmbx(responses,process_units); break; + case RheologyBbarAbsGradientEnum:RheologyBbarAbsGradientx(responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break; case DragCoefficientAbsGradientEnum:DragCoefficientAbsGradientx(responses, elements,nodes, vertices, loads, materials, parameters,process_units,weight_index); break; - case MaterialsRheologyBbarEnum:ElementResponsex(responses, elements,nodes, vertices, loads, materials, parameters,MaterialsRheologyBbarEnum,process_units); break; - case VelEnum:ElementResponsex(responses, elements,nodes, vertices, loads, materials, parameters,VelEnum,process_units); break; + case MaterialsRheologyBbarEnum: this->ElementResponsex(responses,MaterialsRheologyBbarEnum,process_units); break; + case VelEnum: this->ElementResponsex(responses,VelEnum,process_units); break; case FrictionCoefficientEnum:NodalValuex(responses, FrictionCoefficientEnum,elements,nodes, vertices, loads, materials, parameters,process_units); break; default: _error_("response descriptor \"" << EnumToStringx(response_descriptor_enum) << "\" not supported yet!"); break; @@ -897,4 +897,88 @@ }/*}}}*/ +void FemModel::TotalSmbx(IssmDouble* pSmb,bool process_units){/*{{{*/ + + IssmDouble local_smb = 0; + IssmDouble total_smb; + + for(int i=0;ielements->Size();i++){ + Element* element=dynamic_cast(this->elements->GetObjectByOffset(i)); + local_smb+=element->TotalSmb(); + } +#ifdef _HAVE_MPI_ + MPI_Reduce(&local_smb,&total_smb,1,MPI_DOUBLE,MPI_SUM,0,IssmComm::GetComm() ); + MPI_Bcast(&total_smb,1,MPI_DOUBLE,0,IssmComm::GetComm()); +#else + total_smb=local_smb; +#endif + + /*Assign output pointers: */ + *pSmb=total_smb; + +}/*}}}*/ +void FemModel::IceVolumex(IssmDouble* pV,bool process_units){/*{{{*/ + + IssmDouble local_ice_volume = 0; + IssmDouble total_ice_volume; + + for(int i=0;ielements->Size();i++){ + Element* element=dynamic_cast(this->elements->GetObjectByOffset(i)); + local_ice_volume+=element->IceVolume(); + } + #ifdef _HAVE_MPI_ + MPI_Reduce(&local_ice_volume,&total_ice_volume,1,MPI_DOUBLE,MPI_SUM,0,IssmComm::GetComm() ); + MPI_Bcast(&total_ice_volume,1,MPI_DOUBLE,0,IssmComm::GetComm()); + #else + total_ice_volume=local_ice_volume; + #endif + + /*Assign output pointers: */ + *pV=total_ice_volume; + +}/*}}}*/ +void FemModel::ElementResponsex(IssmDouble* presponse,int response_enum,bool process_units){/*{{{*/ + + int found=0; + int sumfound=0; + int cpu_found=-1; + int index; + IssmDouble response; + Element* element=NULL; + + /*retrieve element we are interested in: */ + this->parameters->FindParam(&index,IndexEnum); + int my_rank=IssmComm::GetRank(); + + /*now, go through our elements, and retrieve the one with this id: index: */ + for(int i=0;ielements->Size();i++){ + element=dynamic_cast(this->elements->GetObjectByOffset(i)); + if (element->Id()==index){ + found=1; + cpu_found=my_rank; + break; + } + } + + /*Broadcast whether we found the element: */ +#ifdef _HAVE_MPI_ + MPI_Allreduce ( &found,&sumfound,1,MPI_INT,MPI_SUM,IssmComm::GetComm()); + if(!sumfound)_error_("could not find material with id" << index << " to compute ElementResponse"); +#endif + + /*Ok, we found the element, compute responseocity: */ + if(my_rank==cpu_found){ + element->ElementResponse(&response,response_enum,IuToExtEnum); + } + + /*Broadcast and plug into response: */ +#ifdef _HAVE_MPI_ + MPI_Allreduce ( &cpu_found,&cpu_found,1,MPI_INT,MPI_MAX,IssmComm::GetComm()); + MPI_Bcast(&response,1,MPI_DOUBLE,cpu_found,IssmComm::GetComm()); +#endif + + /*Assign output pointers: */ + *presponse=response; + +}/*}}}*/ #endif #ifdef _HAVE_CONTROL_ Index: /issm/trunk-jpl/src/c/classes/FemModel.h =================================================================== --- /issm/trunk-jpl/src/c/classes/FemModel.h (revision 13874) +++ /issm/trunk-jpl/src/c/classes/FemModel.h (revision 13875) @@ -74,4 +74,7 @@ void MinVyx(IssmDouble* presponse,bool process_units); void MinVzx(IssmDouble* presponse,bool process_units); + void TotalSmbx(IssmDouble* pSmb,bool process_units); + void IceVolumex(IssmDouble* pV,bool process_units); + void ElementResponsex(IssmDouble* presponse,int response_enum,bool process_units); #endif #ifdef _HAVE_DAKOTA_ Index: /issm/trunk-jpl/src/c/modules/modules.h =================================================================== --- /issm/trunk-jpl/src/c/modules/modules.h (revision 13874) +++ /issm/trunk-jpl/src/c/modules/modules.h (revision 13875) @@ -34,6 +34,4 @@ #include "./GroundinglineMigrationx/GroundinglineMigrationx.h" #include "./HoleFillerx/HoleFillerx.h" -#include "./TotalSmbx/TotalSmbx.h" -#include "./IceVolumex/IceVolumex.h" #include "./InputControlUpdatex/InputControlUpdatex.h" #include "./InputConvergencex/InputConvergencex.h" @@ -63,5 +61,4 @@ #include "./Krigingx/Krigingx.h" #include "./Shp2Kmlx/Shp2Kmlx.h" -#include "./ElementResponsex/ElementResponsex.h" #include "./Mergesolutionfromftogx/Mergesolutionfromftogx.h" #include "./MeshPartitionx/MeshPartitionx.h"