diff --git a/autotest/pst_from_tests.py b/autotest/pst_from_tests.py index d7dc0b990..6b44d4d76 100644 --- a/autotest/pst_from_tests.py +++ b/autotest/pst_from_tests.py @@ -509,7 +509,7 @@ def mf6_freyberg_test(): sim.simulation_data.mfpath.set_sim_path(tmp_model_ws) # sim.set_all_data_external() m = sim.get_model("freyberg6") - sim.set_all_data_external() + sim.set_all_data_external(check_data=False) sim.write_simulation() # to by pass the issues with flopy @@ -633,11 +633,22 @@ def mf6_freyberg_test(): datetime=dts[kper]) else: for arr_file in arr_files: + + # these ult bounds are used later in an assert + ult_lb = None + ult_ub = None + if "npf_k_" in arr_file: + ult_ub = 20.0 + ult_lb = 2.0 pf.add_parameters(filenames=arr_file,par_type="grid",par_name_base=arr_file.split('.')[1]+"_gr", pargp=arr_file.split('.')[1]+"_gr",zone_array=ib,upper_bound=ub,lower_bound=lb, - geostruct=gr_gs) + geostruct=gr_gs,ult_ubound=None if ult_ub is None else ult_ub + 1, + ult_lbound=None if ult_lb is None else ult_lb + 1) + # use a slightly lower ult bound here pf.add_parameters(filenames=arr_file, par_type="pilotpoints", par_name_base=arr_file.split('.')[1]+"_pp", - pargp=arr_file.split('.')[1]+"_pp", zone_array=ib,upper_bound=ub,lower_bound=lb,) + pargp=arr_file.split('.')[1]+"_pp", zone_array=ib,upper_bound=ub,lower_bound=lb, + ult_ubound=None if ult_ub is None else ult_ub - 1, + ult_lbound=None if ult_lb is None else ult_lb - 1) # add SP1 spatially constant, but temporally correlated wel flux pars @@ -810,8 +821,6 @@ def mf6_freyberg_test(): print(pst.phi) #assert pst.phi < 1.0e-5, pst.phi - - # check mult files are in pst input files csv = os.path.join(template_ws, "mult2model_info.csv") df = pd.read_csv(csv, index_col=0) @@ -820,6 +829,15 @@ def mf6_freyberg_test(): set(df.loc[df.pp_file.notna()].mlt_file)) assert len(mults_not_linked_to_pst) == 0, print(mults_not_linked_to_pst) + # make sure the appropriate ult bounds have made it thru + df = pd.read_csv(os.path.join(template_ws,"mult2model_info.csv")) + print(df.columns) + df = df.loc[df.model_file.apply(lambda x: "npf_k_" in x),:] + print(df) + print(df.upper_bound) + print(df.lower_bound) + assert np.abs(float(df.upper_bound.min()) - 19.) < 1.0e-6,df.upper_bound.min() + assert np.abs(float(df.lower_bound.max()) - 3.) < 1.0e-6,df.lower_bound.max() def mf6_freyberg_shortnames_test(): import numpy as np @@ -1313,9 +1331,9 @@ def mf6_freyberg_direct_test(): raise Exception("recharge too diff") if __name__ == "__main__": - freyberg_test() + # freyberg_test() # freyberg_prior_build_test() - # mf6_freyberg_test() + mf6_freyberg_test() # mf6_freyberg_shortnames_test() # mf6_freyberg_da_test() # mf6_freyberg_direct_test() diff --git a/pyemu/utils/helpers.py b/pyemu/utils/helpers.py index cdf640ed8..0316bbb32 100644 --- a/pyemu/utils/helpers.py +++ b/pyemu/utils/helpers.py @@ -377,18 +377,18 @@ def calc_observation_ensemble_quantiles(ens, pst, quantiles, subset_obsnames=Non if subset_obsnames is not None: trimnames = subset_obsnames - if len(set(trimnames) - set(obs.index.values)) is not 0: + if len(set(trimnames) - set(obs.index.values)) != 0: raise Exception('the following names in subset_obsnames are not in the ensemble:\n' + ['{}\n'.format(i) for i in (set(trimnames) - set(obs.index.values))]) if subset_obsgroups is not None: - if len((set(subset_obsgroups) - set(pst.obs_groups))) is not 0: + if len((set(subset_obsgroups) - set(pst.obs_groups))) != 0: raise Exception('the following groups in subset_obsgroups are not in pst:\n' + ['{}\n'.format(i) for i in (set(subset_obsgroups) - set(pst.obs_groups))]) trimnames = obs.loc[obs.obgnme.isin(subset_obsgroups)].obsnme.tolist() - if len((set(trimnames) - set(obs.index.values)))is not 0: + if len((set(trimnames) - set(obs.index.values))) != 0: raise Exception('the following names in subset_obsnames are not in the ensemble:\n' + ['{}\n'.format(i) for i in (set(trimnames) - set(obs.index.values))]) # trim the data to subsets (or complete ) diff --git a/pyemu/utils/pst_from.py b/pyemu/utils/pst_from.py index 9cf7a1621..6655f813e 100644 --- a/pyemu/utils/pst_from.py +++ b/pyemu/utils/pst_from.py @@ -122,6 +122,8 @@ def __init__(self, original_d, new_d, longnames=True, self.pst = None self._function_lines_list = [] #each function is itself a list of lines self.direct_org_files = [] + self.ult_ubound_fill = 1.0e+30 + self.ult_lbound_fill = -1.0e+30 @property def parfile_relations(self): @@ -150,25 +152,10 @@ def parfile_relations(self): else pd.Series( {k: v for k, v in [['ubound', x.upper_bound]]}), axis=1) if ubound.nunique(0, False).gt(1).any(): - if ubound.nunique(0, False).gt(2).any(): - # more than one upper bound set - self.logger.lraise( - "different upper bounds requested for same par for {0}" - "".format(name)) - else: - # one set - the rest are None - need to replace None - # with set values - # df with set values - fil = ubound.apply(lambda x: - pd.Series([None]) if x.isna().all() - else x[x.notna()].values).T - self.logger.warn("Upper bound for par passed for some but " - "not all instances, will set NA to " - "passed values\n{}".format(fil)) - # replace Nones in list in Series with passed values - pr.loc[g.index, 'upper_bound'] = g.use_cols.apply( - lambda x: [fil[0].loc['ubound{0}'.format(c)] for c in x] - if x is not None else fil[0].loc['ubound']) + ub_min = ubound.min().fillna(self.ult_ubound_fill).to_dict() + pr.loc[g.index, 'upper_bound'] = g.use_cols.apply( + lambda x: [ub_min['ubound{0}'.format(c)] for c in x] + if x is not None else ub_min["ubound"]) # repeat for lower bounds lbound = g.apply( lambda x: pd.Series( @@ -178,20 +165,10 @@ def parfile_relations(self): else pd.Series( {k: v for k, v in [['lbound', x.lower_bound]]}), axis=1) if lbound.nunique(0, False).gt(1).any(): - if lbound.nunique(0, False).gt(2).any(): - self.logger.lraise( - "different lower bounds requested for same par for {0}" - "".format(name)) - else: - fil = lbound.apply(lambda x: - pd.Series([None]) if x.isna().all() - else x[x.notna()].values).T - self.logger.warn("Lower bound for par passed for some but " - "not all instances, will set NA to " - "passed values\n{}".format(fil)) - pr.loc[g.index, 'lower_bound'] = g.use_cols.apply( - lambda x: [fil[0].loc['lbound{0}'.format(c)] for c in x] - if x is not None else fil[0].loc['lbound']) + lb_max = lbound.max().fillna(self.ult_lbound_fill).to_dict() + pr.loc[g.index, 'lower_bound'] = g.use_cols.apply( + lambda x: [lb_max['lbound{0}'.format(c)] for c in x] + if x is not None else lb_max['lbound']) pr['zero_based'] = self.zero_based return pr @@ -851,6 +828,7 @@ def add_observations(self, filename, insfile=None, use_rows = df.iloc[use_rows].idx_str.unique() # construct ins_file from df ncol = len(use_cols) + obsgp = _check_var_len(obsgp, ncol, fill=True) df_ins = pyemu.pst_utils.csv_to_ins_file( df.set_index('idx_str'), @@ -1048,9 +1026,11 @@ def add_parameters(self, filenames, par_type, zone_array=None, when reading and reapply when writing. Can optionally be `str` in which case `mf_skip` will be treated as a `comment_char`. ult_ubound (`float`): Ultimate upper bound for model input - parameter once all mults are applied - ensure physical model par vals + parameter once all mults are applied - ensure physical model par vals. If not passed, + it is set to 1.0e+30 ult_lbound (`float`): Ultimate lower bound for model input - parameter once all mults are applied + parameter once all mults are applied. If not passed, it is set to + 1.0e-30 for log transform and -1.0e+30 for non-log transform rebuild_pst (`bool`): (Re)Construct PstFrom.pst object after adding new parameters alt_inst_str (`str`): Alternative to default `inst` string in @@ -1070,6 +1050,13 @@ def add_parameters(self, filenames, par_type, zone_array=None, # TODO support passing par_file (i,j)/(x,y) directly where information # is not contained in model parameter file - e.g. no i,j columns + + # this keeps denormal values for creeping into the model input arrays + if ult_ubound is None: + ult_ubound = self.ult_ubound_fill + if ult_lbound is None: + ult_lbound = self.ult_lbound_fill + #some checks for direct parameters par_style = par_style.lower() if par_style not in ["multiplier", "direct"]: @@ -1159,6 +1146,15 @@ def add_parameters(self, filenames, par_type, zone_array=None, "single-element container, or container of " "len use_cols, not '{0}'" "".format(str(par_name_base))) + + # otherewise, things get tripped up in the ensemble/cov stuff + if pargp is not None: + if isinstance(pargp,list): + pargp = [pg.lower() for pg in pargp] + else: + pargp = pargp.lower() + par_name_base = [pnb.lower() for pnb in par_name_base] + if self.longnames: # allow par names to be long... fine for pestpp fmt = "_{0}".format(alt_inst_str) + ":{0}" chk_prefix = "_{0}".format(alt_inst_str) # add `instance` identifier @@ -1169,9 +1165,11 @@ def add_parameters(self, filenames, par_type, zone_array=None, for i in range(len(par_name_base)): par_name_base[i] += fmt.format( self._next_count(par_name_base[i] + chk_prefix)) + # multiplier file name will be taken first par group, if passed # (the same multipliers will apply to all pars passed in this call) # Remove `:` for filenames + par_name_store = par_name_base[0].replace(':', '') # for os filename # Define requisite filenames