Changeset 27226


Ignore:
Timestamp:
08/24/22 13:59:17 (3 years ago)
Author:
Cheng Gong
Message:

CHG: add upper and lower bounds in calving parameterization

Location:
issm/trunk-jpl/src
Files:
6 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp

    r27203 r27226  
    225225                        parameters->AddObject(iomodel->CopyConstantObject("md.calving.xoffset",CalvingXoffsetEnum));
    226226                        parameters->AddObject(iomodel->CopyConstantObject("md.calving.yoffset",CalvingYoffsetEnum));
    227                         parameters->AddObject(iomodel->CopyConstantObject("md.calving.vel_threshold",CalvingVelThresholdEnum));
     227                        parameters->AddObject(iomodel->CopyConstantObject("md.calving.vel_lowerbound",CalvingVelLowerboundEnum));
     228                        parameters->AddObject(iomodel->CopyConstantObject("md.calving.vel_upperbound",CalvingVelUpperboundEnum));
    228229                        break;
    229230                default:

  • issm/trunk-jpl/src/c/classes/Elements/Tria.cpp

    r27203 r27226  
    844844        IssmDouble  arate, rho_ice, rho_water, thickness;
    845845        int                     use_parameter=-1;
    846         IssmDouble  gamma, theta, alpha, xoffset, yoffset, vel_threshold, vrate;
     846        IssmDouble  gamma, theta, alpha, xoffset, yoffset;
     847        IssmDouble  vel_lower, vel_upper, vrate, truncateVrate;
    847848
    848849        /* Get node coordinates and dof list: */
    … …  
    870871        this->FindParam(&xoffset, CalvingXoffsetEnum);
    871872        this->FindParam(&yoffset, CalvingYoffsetEnum);
    872         this->FindParam(&vel_threshold, CalvingVelThresholdEnum);
     873        this->FindParam(&vel_lower, CalvingVelLowerboundEnum);
     874        this->FindParam(&vel_upper, CalvingVelUpperboundEnum);
    873875
    874876        /* Start looping on the number of vertices: */
    … …  
    889891                arate_input->GetInputValue(&arate,gauss);
    890892                vrate = 1.0;
    891                 if (vel < vel_threshold) vrate = vel / vel_threshold;
     893                if (vel < vel_upper) vrate = vel / vel_upper;
    892894
    893895                /*Compute strain rate and viscosity: */
    … …  
    920922                                gamma = yoffset -  0.5*theta*tanh(alpha*(-thickness+xoffset));
    921923                                break;
     924                        case 3:
     925                                /* 3 tanh(normal vel): f(x)=y_{o}-\frac{\theta}{2}\tanh(\alpha(x+x_{o})) */
     926                                _error_("The normalized velocity is not supported yet!");
     927                                gamma = yoffset -  0.5*theta*tanh(alpha*(vel+xoffset));
     928                                break;
     929                        case 4:
     930                                /* 4 tanh(truncated vel): f(x)=y_{o}-\frac{\theta}{2}\tanh(\alpha(x+x_{o})) */
     931                                truncateVrate = (min(vel_upper, max(vel_lower, vel))-vel_lower) / vel_upper;
     932                                gamma = 0.5*theta*(tanh(alpha*xoffset) - tanh(alpha*(truncateVrate+xoffset)));
     933                                break;
    922934                        case -1:
    923935                                /* nothing, just the arate*/
    … …  
    934946
    935947                /*-------------------------------------------*/
    936                 calvingrate[iv] = arate*gamma*vrate;
    937         }
    938 
     948                if (use_parameter < 3) {
     949                        calvingrate[iv] = arate*gamma*vrate;
     950                }
     951                else {
     952                        calvingrate[iv] = arate*gamma;
     953                }
     954        }
    939955        /*Add input*/
    940956        this->AddInput(CalvingCalvingrateEnum,&calvingrate[0],P1DGEnum);

  • issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h

    r27203 r27226  
    116116        CalvingXoffsetEnum,
    117117        CalvingYoffsetEnum,
    118         CalvingVelThresholdEnum,
     118        CalvingVelLowerboundEnum,
     119        CalvingVelUpperboundEnum,
    119120        ConfigurationTypeEnum,
    120121        ConstantsGEnum,

  • issm/trunk-jpl/src/m/classes/calvingparameterization.m

    r27203 r27226  
    1212                xoffset = 0;
    1313                yoffset = 0;
    14                 vel_threshold = 6000;
     14                vel_upperbound = 6000;
     15                vel_lowerbound = 0;
    1516        end
    1617        methods
    … …  
    5253                        % offset in y-axis
    5354                        self.yoffset = 0;
    54                         % velocity threshold to reduce calving rate
    55                         vel_threshold = 6000; % m/a
     55                        % velocity thresholds to reduce calving rate
     56                        vel_upperbound = 6000; % m/a
     57                        vel_lowerbound = 0; % m/a
    5658                end % }}}
    5759                function md = checkconsistency(self,md,solution,analyses) % {{{
    … …  
    6062
    6163                        md = checkfield(md,'fieldname','calving.min_thickness','>=',0,'NaN',1,'Inf',1,'numel',1);
    62                         md = checkfield(md,'fieldname','calving.use_param','values',[-1, 0, 1, 2]);
     64                        md = checkfield(md,'fieldname','calving.use_param','values',[-1, 0, 1, 2, 3, 4]);
    6365                        md = checkfield(md,'fieldname','calving.theta','NaN',1,'Inf',1,'numel',1);
    6466                        md = checkfield(md,'fieldname','calving.alpha','NaN',1,'Inf',1,'numel',1);
    6567                        md = checkfield(md,'fieldname','calving.xoffset','NaN',1,'Inf',1,'numel',1);
    6668                        md = checkfield(md,'fieldname','calving.yoffset','NaN',1,'Inf',1,'numel',1);
    67                         md = checkfield(md,'fieldname','calving.vel_threshold','NaN',1,'Inf',1,'numel',1);
     69                        md = checkfield(md,'fieldname','calving.vel_lowerbound','NaN',1,'Inf',1,'numel',1);
     70                        md = checkfield(md,'fieldname','calving.vel_upperbound','NaN',1,'Inf',1,'numel',1);
    6871                end % }}}
    6972                function disp(self) % {{{
    7073                        disp(sprintf('   Calving test parameters:'));
    7174                        fielddisplay(self,'min_thickness','minimum thickness below which no ice is allowed [m]');
    72                         fielddisplay(self,'use_param','-1 - just use frontal ablation rate, 0 - f(x) = y_{o} + \alpha (x+x_{o}), 1 - f(x)=y_{o}-\frac{\theta}{2}\tanh(\alpha(x+x_{o})), 2 - tanh(thickness)');
     75                        fielddisplay(self,'use_param','-1 - just use frontal ablation rate, 0 - f(x) = y_{o} + \alpha (x+x_{o}), 1 - f(x)=y_{o}-\frac{\theta}{2}\tanh(\alpha(x+x_{o})), 2 - tanh(thickness), 3 - tanh(normalized vel), 4 - tanh(truncated vel)');
    7376                        fielddisplay(self,'theta','the amplifier');
    7477                        fielddisplay(self,'alpha','the slope');
    7578                        fielddisplay(self,'xoffset','offset in x-axis');
    7679                        fielddisplay(self,'yoffset','offset in y-axis');
    77                         fielddisplay(self,'vel_threshold','threshold of ice velocity to reduce the calving rate');
     80                        fielddisplay(self,'vel_lowerbound','lowerbound of ice velocity to reduce the calving rate [m/a]');
     81                        fielddisplay(self,'vel_upperbound','upperbound of ice velocity to reduce the calving rate [m/a]');
    7882                end % }}}
    7983                function marshall(self,prefix,md,fid) % {{{
    … …  
    8690                        WriteData(fid,prefix,'object',self,'fieldname','xoffset','format','Double');
    8791                        WriteData(fid,prefix,'object',self,'fieldname','yoffset','format','Double');
    88                         WriteData(fid,prefix,'object',self,'fieldname','vel_threshold','format','Double','scale', 1./yts);
     92                        WriteData(fid,prefix,'object',self,'fieldname','vel_lowerbound','format','Double','scale', 1./yts);
     93                        WriteData(fid,prefix,'object',self,'fieldname','vel_upperbound','format','Double','scale', 1./yts);
    8994                end % }}}
    9095        end

  • issm/trunk-jpl/src/m/classes/initialization.m

    r27021 r27226  
    150150                        yts=md.constants.yts;
    151151
    152                         WriteData(fid,prefix,'object',self,'fieldname','vx','format','DoubleMat','mattype',1,'scale',1./yts);
    153                         WriteData(fid,prefix,'object',self,'fieldname','vy','format','DoubleMat','mattype',1,'scale',1./yts);
     152                        WriteData(fid,prefix,'object',self,'fieldname','vx','format','DoubleMat','mattype',1,'scale',1./yts,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);
     153                        WriteData(fid,prefix,'object',self,'fieldname','vy','format','DoubleMat','mattype',1,'scale',1./yts,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);
    154154                        WriteData(fid,prefix,'object',self,'fieldname','vz','format','DoubleMat','mattype',1,'scale',1./yts);
    155155                        WriteData(fid,prefix,'object',self,'fieldname','pressure','format','DoubleMat','mattype',1);

  • issm/trunk-jpl/src/m/parameterization/setflowequation.m

    r27031 r27226  
    5555%check that each element has at least one flag
    5656if any(SIAflag+SSAflag+HOflag+L1L2flag+MOLHOflag+FSflag==0),
    57         error('elements type not assigned, supported models are ''SIA'',''SSA'',''HO'' and ''FS''')
     57        error('elements type not assigned, supported models are ''SIA'',''SSA'',''HO'',''MOLHO'' and ''FS''')
    5858end
    5959
Note: See TracChangeset for help on using the changeset viewer.