Changeset 23077

Timestamp:

08/07/18 21:46:18 (7 years ago)

Author:

seroussi

Message:

CHG: make sure base is +9999. in places where there is no more ice

File:

• issm/trunk-jpl/src/c/cores/transient_core.cpp (modified) (8 diffs)

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

@@ -93 +93 @@
                 IssmDouble *oceangridx;
                 IssmDouble *oceangridy;
-                IssmDouble *base_oceangrid;
+                IssmDouble *icebase_oceangrid;
                 IssmDouble* x_ice = NULL;
                 IssmDouble* y_ice = NULL;
@@ -153 +153 @@
                 odouble->size=NULL;
                 options->AddOption(odouble);
-                InterpFromMeshToMesh2dx(&base_oceangrid,index_ice,lon_ice,lat_ice,ngrids_ice,nels_ice,
+                InterpFromMeshToMesh2dx(&icebase_oceangrid,index_ice,lon_ice,lat_ice,ngrids_ice,nels_ice,
                                                 icebase,ngrids_ice,1,
                                                 oceangridx,oceangridy,ngrids_ocean,options);
@@ -159 +159 @@
 
                 if(my_rank==0){
-                        ISSM_MPI_Send(base_oceangrid,ngrids_ocean,ISSM_MPI_DOUBLE,0,10001008,tomitgcmcomm);
+                        ISSM_MPI_Send(icebase_oceangrid,ngrids_ocean,ISSM_MPI_DOUBLE,0,10001008,tomitgcmcomm);
                 }
 
@@ -171 +171 @@
                 xDelete<IssmDouble>(z_ice);
                 xDelete<IssmDouble>(icebase);
-                xDelete<IssmDouble>(base_oceangrid);
+                xDelete<IssmDouble>(icebase_oceangrid);
                 xDelete<IssmDouble>(oceangridx);
                 xDelete<IssmDouble>(oceangridy);
@@ -237 +237 @@
                         IssmDouble oceantime;
                         IssmDouble rho_ice;
-                        IssmDouble *oceanmelt;
-                        IssmDouble *base_oceangrid;
-                        IssmDouble *oceangridx;
-                        IssmDouble *oceangridy;
-                        IssmDouble* x_ice = NULL;
-                        IssmDouble* y_ice = NULL;
-                        IssmDouble* lat_ice = NULL;
-                        IssmDouble* lon_ice = NULL;
-                        IssmDouble* z_ice = NULL;
-                        IssmDouble* icebase= NULL;
-                        IssmDouble* melt_mesh = NULL;
-                        int*        index_ice= NULL;
-                        int*        index_ocean = NULL;
+                        IssmDouble *oceanmelt         = NULL;
+                        IssmDouble *icebase_oceangrid = NULL;
+                        IssmDouble *icemask_oceangrid = NULL;
+                        IssmDouble *oceangridx        = NULL;
+                        IssmDouble *oceangridy        = NULL;
+                        IssmDouble *x_ice             = NULL;
+                        IssmDouble *y_ice             = NULL;
+                        IssmDouble *lat_ice           = NULL;
+                        IssmDouble *lon_ice           = NULL;
+                        IssmDouble *z_ice             = NULL;
+                        IssmDouble *icebase           = NULL;
+                        IssmDouble *icemask           = NULL;
+                        IssmDouble *melt_mesh         = NULL;
+                        int        *index_ice         = NULL;
+                        int        *index_ocean       = NULL;
                         int         ngrids_ice=femmodel->vertices->NumberOfVertices();
                         int         nels_ice=femmodel->elements->NumberOfElements();
@@ -262 +264 @@
                         femmodel->vertices->LatLonList(&lat_ice,&lon_ice);
 
-                        /*Interpolate ice base onto ocean grid*/
+                        /*Interpolate ice base and mask onto ocean grid*/
                         GetVectorFromInputsx(&icebase,femmodel,BaseEnum,VertexSIdEnum);
-                        InterpFromMeshToMesh2dx(&base_oceangrid,index_ice,lon_ice,lat_ice,ngrids_ice,nels_ice,
-                                                icebase,ngrids_ice,1,
-                                                oceangridx,oceangridy,ngrids_ocean,NULL);
+                        InterpFromMeshToMesh2dx(&icebase_oceangrid,index_ice,lon_ice,lat_ice,ngrids_ice,nels_ice,
+                                                icebase,ngrids_ice,1,oceangridx,oceangridy,ngrids_ocean,NULL);
+                        xDelete<IssmDouble>(icebase);
+
+                        GetVectorFromInputsx(&icemask,femmodel,BaseEnum,VertexSIdEnum);
+                        InterpFromMeshToMesh2dx(&icemask_oceangrid,index_ice,lon_ice,lat_ice,ngrids_ice,nels_ice,
+                                                icemask,ngrids_ice,1,oceangridx,oceangridy,ngrids_ocean,NULL);
+                        xDelete<IssmDouble>(icemask);
+
+                        /*Put +9999 for places where there is no ice!*/
+                        for(int i=0;i<ngrids_ocean;i++) if(icemask_oceangrid[i]>0.) icebase_oceangrid[i]=+9999.;
+                        xDelete<IssmDouble>(icemask_oceangrid);
 
                         /*Send and receive data*/
@@ -275 +286 @@
                                 oceanmelt = xNew<IssmDouble>(ngrids_ocean);
                                 ISSM_MPI_Recv(oceanmelt,ngrids_ocean,ISSM_MPI_DOUBLE,0,10001007,tomitgcmcomm,&status);
-                                ISSM_MPI_Send(base_oceangrid,ngrids_ocean,ISSM_MPI_DOUBLE,0,10001008,tomitgcmcomm);
+                                ISSM_MPI_Send(icebase_oceangrid,ngrids_ocean,ISSM_MPI_DOUBLE,0,10001008,tomitgcmcomm);
                         }
                         ISSM_MPI_Bcast(&oceantime,1,ISSM_MPI_DOUBLE,0,IssmComm::GetComm());
@@ -301 +312 @@
                         xDelete<IssmDouble>(oceangridy);
                         xDelete<IssmDouble>(oceanmelt);
-                        xDelete<IssmDouble>(icebase);
-                        xDelete<IssmDouble>(base_oceangrid);
+                        xDelete<IssmDouble>(icebase_oceangrid);
                 }
                 #else