CalvingAtThickness.cc (3597B)
1 /* Copyright (C) 2013, 2014, 2015, 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 "CalvingAtThickness.hh" 21 22 #include "pism/util/Mask.hh" 23 #include "pism/util/IceGrid.hh" 24 #include "pism/util/pism_utilities.hh" 25 26 namespace pism { 27 28 //! @brief Calving and iceberg removal code. 29 namespace calving { 30 31 CalvingAtThickness::CalvingAtThickness(IceGrid::ConstPtr g) 32 : Component(g), 33 m_calving_threshold(m_grid, "thickness_calving_threshold", WITHOUT_GHOSTS), 34 m_old_mask(m_grid, "old_mask", WITH_GHOSTS, 1) { 35 36 m_calving_threshold.set_attrs("diagnostic", 37 "threshold used by the 'calving at threshold' calving method", 38 "m", "m", 39 "", 0); // no standard name 40 m_calving_threshold.set_time_independent(true); 41 } 42 43 CalvingAtThickness::~CalvingAtThickness() { 44 // empty 45 } 46 47 48 void CalvingAtThickness::init() { 49 50 m_log->message(2, "* Initializing the 'calving at a threshold thickness' mechanism...\n"); 51 52 std::string threshold_file = m_config->get_string("calving.thickness_calving.threshold_file"); 53 54 double calving_threshold = m_config->get_number("calving.thickness_calving.threshold"); 55 56 m_calving_threshold.set(calving_threshold); 57 58 if (not threshold_file.empty()) { 59 m_log->message(2, 60 " Reading thickness calving threshold from file '%s'...\n", 61 threshold_file.c_str()); 62 63 m_calving_threshold.regrid(threshold_file, CRITICAL); 64 } else { 65 m_log->message(2, 66 " Thickness threshold: %3.3f meters.\n", calving_threshold); 67 } 68 } 69 70 /** 71 * Updates ice cover mask and the ice thickness according to the 72 * calving rule removing ice at the shelf front that is thinner than a 73 * given threshold. 74 * 75 * @param[in,out] pism_mask ice cover mask 76 * @param[in,out] ice_thickness ice thickness 77 * 78 * @return 0 on success 79 */ 80 void CalvingAtThickness::update(IceModelVec2CellType &pism_mask, 81 IceModelVec2S &ice_thickness) { 82 83 // this call fills ghosts of m_old_mask 84 m_old_mask.copy_from(pism_mask); 85 86 IceModelVec::AccessList list{&pism_mask, &ice_thickness, &m_old_mask, &m_calving_threshold}; 87 for (Points p(*m_grid); p; p.next()) { 88 const int i = p.i(), j = p.j(); 89 90 if (m_old_mask.floating_ice(i, j) && 91 m_old_mask.next_to_ice_free_ocean(i, j) && 92 ice_thickness(i, j) < m_calving_threshold(i, j)) { 93 ice_thickness(i, j) = 0.0; 94 pism_mask(i, j) = MASK_ICE_FREE_OCEAN; 95 } 96 } 97 98 pism_mask.update_ghosts(); 99 ice_thickness.update_ghosts(); 100 } 101 102 const IceModelVec2S& CalvingAtThickness::threshold() const { 103 return m_calving_threshold; 104 } 105 106 DiagnosticList CalvingAtThickness::diagnostics_impl() const { 107 return {{"thickness_calving_threshold", Diagnostic::wrap(m_calving_threshold)}}; 108 } 109 110 } // end of namespace calving 111 } // end of namespace pism