Changeset 23714


Ignore:
Timestamp:
02/11/19 20:14:02 (6 years ago)
Author:
seroussi
Message:

CHG: cleaning up ocean coupling

Location:
issm/trunk-jpl/src/c
Files:
3 added
3 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk-jpl/src/c/Makefile.am

    r23652 r23714  
    520520endif
    521521#}}}
     522#Oceansources  {{{
     523if OCEAN
     524issm_sources +=  ./modules/OceanExchangeDatax/OceanExchangeDatax.cpp
     525endif
     526#}}}
    522527#Slr sources  {{{
    523528if SEALEVELRISE

  • issm/trunk-jpl/src/c/cores/transient_core.cpp

    r23665 r23714  
    7979        #endif
    8080
    81         if(isoceancoupling){ /*{{{*/
    82                 #ifndef _HAVE_AD_
    83                 if(VerboseSolution()) _printf0_("   ocean coupling: initialization \n");
    84                 int my_rank;
    85                 ISSM_MPI_Comm tomitgcmcomm;
    86                 ISSM_MPI_Status status;
    87 
    88                 my_rank=IssmComm::GetRank();
    89                 GenericParam<ISSM_MPI_Comm>* parcom = dynamic_cast<GenericParam<ISSM_MPI_Comm>*>(femmodel->parameters->FindParamObject(ToMITgcmCommEnum));
    90                 if(!parcom)_error_("TransferForcing error message: could not find ToMITgcmCommEnum communicator");
    91                 tomitgcmcomm=parcom->GetParameterValue();
    92 
    93                 int oceangridnxsize,oceangridnysize,ngrids_ocean,nels_ocean;
    94                 IssmDouble  oceantime,coupling_time;
    95                 IssmDouble *oceangridx;
    96                 IssmDouble *oceangridy;
    97                 IssmDouble *icebase_oceangrid = NULL;
    98                 IssmDouble *icemask_oceangrid = NULL;
    99                 IssmDouble* x_ice             = NULL;
    100                 IssmDouble* y_ice             = NULL;
    101                 IssmDouble* lat_ice           = NULL;
    102                 IssmDouble* lon_ice           = NULL;
    103                 IssmDouble* icebase           = NULL;
    104                 IssmDouble* icemask           = NULL;
    105                 int*        index_ice         = NULL;
    106                 int*        index_ocean       = NULL;
    107                 int         ngrids_ice=femmodel->vertices->NumberOfVertices();
    108                 int         nels_ice=femmodel->elements->NumberOfElements();
    109 
    110                 /*Recover fixed parameters and store them*/
    111                 femmodel->parameters->FindParam(&coupling_time,TimesteppingCouplingTimeEnum);
    112                 if(my_rank==0){
    113                         ISSM_MPI_Send(&coupling_time,1,ISSM_MPI_DOUBLE,0,10001000,tomitgcmcomm);
    114                         ISSM_MPI_Recv(&oceangridnxsize,1,ISSM_MPI_INT,0,10001003,tomitgcmcomm,&status);
    115                         ISSM_MPI_Recv(&oceangridnysize,1,ISSM_MPI_INT,0,10001004,tomitgcmcomm,&status);
    116                 }
    117                 ngrids_ocean=oceangridnxsize*oceangridnysize;
    118                 ISSM_MPI_Bcast(&oceangridnxsize,1,ISSM_MPI_INT,0,IssmComm::GetComm());
    119                 ISSM_MPI_Bcast(&oceangridnysize,1,ISSM_MPI_INT,0,IssmComm::GetComm());
    120                 ISSM_MPI_Bcast(&ngrids_ocean,1,ISSM_MPI_INT,0,IssmComm::GetComm());
    121                 ISSM_MPI_Bcast(&oceantime,1,ISSM_MPI_DOUBLE,0,IssmComm::GetComm());
    122                 femmodel->parameters->SetParam(oceangridnxsize,OceanGridNxEnum);
    123                 femmodel->parameters->SetParam(oceangridnysize,OceanGridNyEnum);
    124                 if(my_rank==0){
    125                         oceangridx = xNew<IssmDouble>(ngrids_ocean);
    126                         ISSM_MPI_Recv(oceangridx,ngrids_ocean,ISSM_MPI_DOUBLE,0,10001005,tomitgcmcomm,&status);
    127                         oceangridy = xNew<IssmDouble>(ngrids_ocean);
    128                         ISSM_MPI_Recv(oceangridy,ngrids_ocean,ISSM_MPI_DOUBLE,0,10001006,tomitgcmcomm,&status);
    129 
    130                         /*Exchange varying parameters for the initialization*/
    131                         ISSM_MPI_Send(&time,1,ISSM_MPI_DOUBLE,0,10001001,tomitgcmcomm);
    132                         ISSM_MPI_Recv(&oceantime,1,ISSM_MPI_DOUBLE,0,10001002,tomitgcmcomm,&status);
    133                 }
    134                 if(my_rank!=0){
    135                         oceangridx=xNew<IssmDouble>(ngrids_ocean);
    136                         oceangridy=xNew<IssmDouble>(ngrids_ocean);
    137                 }
    138                 ISSM_MPI_Bcast(oceangridx,ngrids_ocean,ISSM_MPI_DOUBLE,0,IssmComm::GetComm());
    139                 ISSM_MPI_Bcast(oceangridy,ngrids_ocean,ISSM_MPI_DOUBLE,0,IssmComm::GetComm());
    140                 femmodel->parameters->SetParam(oceangridx,ngrids_ocean,OceanGridXEnum);
    141                 femmodel->parameters->SetParam(oceangridy,ngrids_ocean,OceanGridYEnum);
    142 
    143                 /*Interpolate ice base and mask onto ocean grid*/
    144                 femmodel->GetMesh(femmodel->vertices,femmodel->elements,&x_ice,&y_ice,&index_ice);
    145                 BamgTriangulatex(&index_ocean,&nels_ocean,oceangridx,oceangridy,ngrids_ocean);
    146                 femmodel->vertices->LatLonList(&lat_ice,&lon_ice);
    147                 GetVectorFromInputsx(&icebase,femmodel,BaseEnum,VertexSIdEnum);
    148                 Options* options = new Options();
    149                 GenericOption<double> *odouble = new GenericOption<double>();
    150                 const char* name = "default";
    151                 odouble->name =xNew<char>(strlen(name)+1);
    152                 memcpy(odouble->name,name,(strlen(name)+1)*sizeof(char));
    153                 odouble->value=+9999.;
    154                 odouble->size[0]=1;
    155                 odouble->size[1]=1;
    156                 options->AddOption(odouble);
    157                 InterpFromMeshToMesh2dx(&icebase_oceangrid,index_ice,lon_ice,lat_ice,ngrids_ice,nels_ice,
    158                                                 icebase,ngrids_ice,1,
    159                                                 oceangridx,oceangridy,ngrids_ocean,options);
    160                 delete options;
    161                 xDelete<IssmDouble>(icebase);
    162 
    163                 GetVectorFromInputsx(&icemask,femmodel,MaskIceLevelsetEnum,VertexSIdEnum);
    164                 Options* options2 = new Options();
    165                 GenericOption<double> *odouble2 = new GenericOption<double>();
    166                 const char* name2 = "default";
    167                 odouble2->name =xNew<char>(strlen(name2)+1);
    168                 memcpy(odouble2->name,name2,(strlen(name2)+1)*sizeof(char));
    169                 odouble2->value=+1.;
    170                 odouble2->size[0]=1;
    171                 odouble2->size[1]=1;
    172                 options2->AddOption(odouble2);
    173                 InterpFromMeshToMesh2dx(&icemask_oceangrid,index_ice,lon_ice,lat_ice,ngrids_ice,nels_ice,
    174                                         icemask,ngrids_ice,1,oceangridx,oceangridy,ngrids_ocean,options2);
    175                 delete options2;
    176                 xDelete<IssmDouble>(icemask);
    177 
    178                 /*Put +9999 for places where there is no ice!*/
    179                 for(int i=0;i<ngrids_ocean;i++) if(icemask_oceangrid[i]>0.) icebase_oceangrid[i]=+9999.;
    180                 xDelete<IssmDouble>(icemask_oceangrid);
    181 
    182                 if(my_rank==0){
    183                         ISSM_MPI_Send(icebase_oceangrid,ngrids_ocean,ISSM_MPI_DOUBLE,0,10001008,tomitgcmcomm);
    184                 }
    185 
    186                 /*Delete*/
    187                 xDelete<int>(index_ice);
    188                 xDelete<int>(index_ocean);
    189                 xDelete<IssmDouble>(lat_ice);
    190                 xDelete<IssmDouble>(lon_ice);
    191                 xDelete<IssmDouble>(x_ice);
    192                 xDelete<IssmDouble>(y_ice);
    193                 xDelete<IssmDouble>(icebase_oceangrid);
    194                 xDelete<IssmDouble>(oceangridx);
    195                 xDelete<IssmDouble>(oceangridy);
    196         #else
    197         _error_("not supported");
    198         #endif
    199         }
    200         /*}}}*/
    201 
    202                 IssmDouble  output_value;
    203                 int         num_dependents;
    204                 IssmPDouble *dependents;
    205                 DataSet*    dependent_objects=NULL;
    206                 IssmDouble  J=0.;
     81        if(isoceancoupling) OceanExchangeDatax(femmodel);
     82
     83        IssmDouble  output_value;
     84        int         num_dependents;
     85        IssmPDouble *dependents;
     86        DataSet*    dependent_objects=NULL;
     87        IssmDouble  J=0.;
    20788
    20889        if(iscontrol && isautodiff){

  • issm/trunk-jpl/src/c/modules/modules.h

    r23652 r23714  
    100100#include "./UpdateDynamicConstraintsx/UpdateDynamicConstraintsx.h"
    101101#include "./VertexCoordinatesx/VertexCoordinatesx.h"
     102
     103#ifdef _HAVE_OCEAN_
     104#include "./OceanExchangeDatax/OceanExchangeDatax.h"
    102105#endif
     106
     107#endif
Note: See TracChangeset for help on using the changeset viewer.