URI:
       tInitialization.cc - pism - [fork] customized build of PISM, the parallel ice sheet model (tillflux branch)
  HTML git clone git://src.adamsgaard.dk/pism
   DIR Log
   DIR Files
   DIR Refs
   DIR LICENSE
       ---
       tInitialization.cc (6802B)
       ---
            1 /* Copyright (C) 2016, 2017, 2018, 2019 PISM Authors
            2  *
            3  * This file is part of PISM.
            4  *
            5  * PISM is free software; you can redistribute it and/or modify it under the
            6  * terms of the GNU General Public License as published by the Free Software
            7  * Foundation; either version 3 of the License, or (at your option) any later
            8  * version.
            9  *
           10  * PISM is distributed in the hope that it will be useful, but WITHOUT ANY
           11  * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
           12  * FOR A PARTICULAR PURPOSE.  See the GNU General Public License for more
           13  * details.
           14  *
           15  * You should have received a copy of the GNU General Public License
           16  * along with PISM; if not, write to the Free Software
           17  * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301  USA
           18  */
           19 
           20 #include "Initialization.hh"
           21 #include "pism/util/error_handling.hh"
           22 #include "pism/util/pism_utilities.hh"
           23 #include "pism/util/io/File.hh"
           24 #include "pism/coupler/util/init_step.hh"
           25 
           26 namespace pism {
           27 namespace surface {
           28 
           29 InitializationHelper::InitializationHelper(IceGrid::ConstPtr grid, std::shared_ptr<SurfaceModel> input)
           30   : SurfaceModel(grid, input) {
           31 
           32   if (not input) {
           33     throw RuntimeError(PISM_ERROR_LOCATION, "pism::surface::InitializationHelper got a NULL input model");
           34   }
           35 
           36   // allocate storage
           37   {
           38     m_mass_flux.create(m_grid, "effective_climatic_mass_balance", WITHOUT_GHOSTS);
           39     m_mass_flux.set_attrs("model_state",
           40                           "surface mass balance (accumulation/ablation) rate, as seen by the ice dynamics code (used for restarting)",
           41                           "kg m-2 s-1", "kg m-2 year-1", "", 0);
           42     m_mass_flux.set_time_independent(false);
           43 
           44     m_temperature.create(m_grid, "effective_ice_surface_temp", WITHOUT_GHOSTS);
           45     m_temperature.set_attrs("model_state",
           46                             "temperature of the ice at the ice surface but below firn processes, as seen by the ice dynamics code (used for restarting)",
           47                             "Kelvin", "Kelvin", "", 0);
           48     m_temperature.set_time_independent(false);
           49 
           50     m_liquid_water_fraction = allocate_liquid_water_fraction(grid);
           51     m_liquid_water_fraction->metadata().set_name("effective_ice_surface_liquid_water_fraction");
           52     m_liquid_water_fraction->set_attrs("model_state",
           53                                        "liquid water fraction of the ice at the top surface, as seen by the ice dynamics code (used for restarting)",
           54                                        "1", "1", "", 0);
           55     m_liquid_water_fraction->set_time_independent(false);
           56 
           57     m_layer_mass = allocate_layer_mass(grid);
           58     m_layer_mass->metadata().set_name("effective_surface_layer_mass");
           59     m_layer_mass->set_attrs("model_state",
           60                             "mass held in the surface layer, as seen by the ice dynamics code (used for restarting)",
           61                             "kg", "kg",
           62                             "", 0);
           63     m_layer_mass->set_time_independent(false);
           64 
           65     m_layer_thickness = allocate_layer_thickness(grid);
           66     m_layer_thickness->metadata().set_name("effective_surface_layer_thickness");
           67     m_layer_thickness->set_attrs("model_state",
           68                                  "thickness of the surface layer, as seen by the ice dynamics code (used for restarting)",
           69                                  "meters", "meters", "", 0);
           70     m_layer_thickness->set_time_independent(false);
           71   }
           72 
           73   m_accumulation = allocate_accumulation(grid);
           74   m_accumulation->set_name("effective_" + m_accumulation->get_name());
           75 
           76   m_melt = allocate_melt(grid);
           77   m_melt->set_name("effective_" + m_melt->get_name());
           78 
           79   m_runoff = allocate_runoff(grid);
           80   m_runoff->set_name("effective_" + m_runoff->get_name());
           81 
           82   // collect pointers
           83   m_variables = {&m_mass_flux,
           84                  &m_temperature,
           85                  m_liquid_water_fraction.get(),
           86                  m_layer_mass.get(),
           87                  m_layer_thickness.get(),
           88                  m_accumulation.get(),
           89                  m_melt.get(),
           90                  m_runoff.get()};
           91 }
           92 
           93 void InitializationHelper::init_impl(const Geometry &geometry) {
           94   m_input_model->init(geometry);
           95 
           96   InputOptions opts = process_input_options(m_grid->com, m_config);
           97 
           98   if (opts.type == INIT_RESTART) {
           99     m_log->message(2, "* Reading effective surface model outputs from '%s' for re-starting...\n",
          100                    opts.filename.c_str());
          101 
          102     File file(m_grid->com, opts.filename, PISM_GUESS, PISM_READONLY);
          103     const unsigned int last_record = file.nrecords() - 1;
          104     for (auto v : m_variables) {
          105       v->read(file, last_record);
          106     }
          107   } else {
          108     m_log->message(2, "* Performing a 'fake' surface model time-step for bootstrapping...\n");
          109 
          110     init_step(this, geometry, *m_grid->ctx()->time());
          111   }
          112 
          113   // Support regridding. This is needed to ensure that initialization using "-i" is equivalent to
          114   // "-i ... -bootstrap -regrid_file ..."
          115   for (auto v : m_variables) {
          116     regrid("surface model initialization helper", *v, REGRID_WITHOUT_REGRID_VARS);
          117   }
          118 }
          119 
          120 void InitializationHelper::update_impl(const Geometry &geometry, double t, double dt) {
          121   m_input_model->update(geometry, t, dt);
          122 
          123   // store outputs of the input model
          124   m_mass_flux.copy_from(m_input_model->mass_flux());
          125   m_temperature.copy_from(m_input_model->temperature());
          126   m_liquid_water_fraction->copy_from(m_input_model->liquid_water_fraction());
          127   m_layer_mass->copy_from(m_input_model->layer_mass());
          128   m_layer_thickness->copy_from(m_input_model->layer_thickness());
          129   m_accumulation->copy_from(m_input_model->accumulation());
          130   m_melt->copy_from(m_input_model->melt());
          131   m_runoff->copy_from(m_input_model->runoff());
          132 }
          133 
          134 const IceModelVec2S &InitializationHelper::layer_thickness_impl() const {
          135   return *m_layer_thickness;
          136 }
          137 
          138 const IceModelVec2S &InitializationHelper::mass_flux_impl() const {
          139   return m_mass_flux;
          140 }
          141 
          142 const IceModelVec2S &InitializationHelper::temperature_impl() const {
          143   return m_temperature;
          144 }
          145 
          146 const IceModelVec2S &InitializationHelper::liquid_water_fraction_impl() const {
          147   return *m_liquid_water_fraction;
          148 }
          149 
          150 const IceModelVec2S &InitializationHelper::layer_mass_impl() const {
          151   return *m_layer_mass;
          152 }
          153 
          154 const IceModelVec2S &InitializationHelper::accumulation_impl() const {
          155   return *m_accumulation;
          156 }
          157 
          158 const IceModelVec2S &InitializationHelper::melt_impl() const {
          159   return *m_melt;
          160 }
          161 
          162 const IceModelVec2S &InitializationHelper::runoff_impl() const {
          163   return *m_runoff;
          164 }
          165 
          166 void InitializationHelper::define_model_state_impl(const File &output) const {
          167   for (auto v : m_variables) {
          168     v->define(output);
          169   }
          170   m_input_model->define_model_state(output);
          171 }
          172 
          173 void InitializationHelper::write_model_state_impl(const File &output) const {
          174   for (auto v : m_variables) {
          175     v->write(output);
          176   }
          177   m_input_model->write_model_state(output);
          178 }
          179 
          180 
          181 
          182 } // end of namespace surface
          183 } // end of namespace pism