From 3b7a49f81e8f3ac70f139e95b05789aa0376dba0 Mon Sep 17 00:00:00 2001 From: Melissa Rasmussen Date: Thu, 5 Feb 2026 20:52:23 -0500 Subject: [PATCH 01/10] Set up old file format reading in --- read_model.H | 47 +++++++++++++++++++++++++++++++---------------- 1 file changed, 31 insertions(+), 16 deletions(-) diff --git a/read_model.H b/read_model.H index b1a35d2..8c1a2e7 100644 --- a/read_model.H +++ b/read_model.H @@ -157,29 +157,44 @@ read_file(const std::string filename, initial_model_t& initial_model) { std::string line; - // first the header line -- this tells us the number of points + // first the header line -- this will determine the format. If it + // tells us the number of points (npts), then we are the old style. getline(mf, line); - std::string npts_string = line.substr(line.find("=")+1, line.length()); - initial_model.npts = std::stoi(npts_string); - if (initial_model.npts > NPTS_MODEL) { - amrex::Error("Error: model has more than NPTS_MODEL points,"); - } - - // next line tells use the number of variables + if (line.find("npts") != std::string::npos) { + // old format + + std::string npts_string = line.substr(line.find("=")+1, line.length()); + initial_model.npts = std::stoi(npts_string); - getline(mf, line); - std::string num_vars_string = line.substr(line.find("=")+1, line.length()); - int nvars_model_file = std::stoi(num_vars_string); + if (initial_model.npts > NPTS_MODEL) { + amrex::Error("Error: model has more than NPTS_MODEL points,"); + } - // now read in the names of the variables + // next line tells use the number of variables - std::vector varnames_stored; - for (int n = 0; n < nvars_model_file; n++) { getline(mf, line); - std::string var_string = line.substr(line.find("#")+1, line.length()); - varnames_stored.push_back(model_string::ltrim(model_string::rtrim(var_string))); + + if (line.find("num of variables") == std::string::npos) { + amrex::Error("Error: missing the number of variables header"); + } + + std::string num_vars_string = line.substr(line.find("=")+1, line.length()); + int nvars_model_file = std::stoi(num_vars_string); + + // now read in the names of the variables + + std::vector varnames_stored; + for (int n = 0; n < nvars_model_file; n++) { + getline(mf, line); + std::string var_string = line.substr(line.find("#")+1, line.length()); + varnames_stored.push_back(model_string::ltrim(model_string::rtrim(var_string))); + } + } else { + // new format + + } From 3ef065757c4dd8ddae1629df68a6821de309db3b Mon Sep 17 00:00:00 2001 From: Melissa Rasmussen Date: Thu, 5 Feb 2026 21:02:58 -0500 Subject: [PATCH 02/10] Set up new file format, and examples of old and new --- read_model.H | 47 ++++++++++++++++++++++++++++++++++++++++++----- 1 file changed, 42 insertions(+), 5 deletions(-) diff --git a/read_model.H b/read_model.H index 8c1a2e7..9eb37ea 100644 --- a/read_model.H +++ b/read_model.H @@ -162,8 +162,23 @@ read_file(const std::string filename, initial_model_t& initial_model) { getline(mf, line); + // variables: + // initial_model.npts + int nvars_model_file; + std::vector varnames_stored; + if (line.find("npts") != std::string::npos) { - // old format + // old format. example: + """ + # npts = 3 + # number of variables = 2 + # radius + # density + # temperature + 0.1 3 4.0 + 0.2 4 3.9 + 0.3 4 3.7 + """ std::string npts_string = line.substr(line.find("=")+1, line.length()); initial_model.npts = std::stoi(npts_string); @@ -181,20 +196,42 @@ read_file(const std::string filename, initial_model_t& initial_model) { } std::string num_vars_string = line.substr(line.find("=")+1, line.length()); - int nvars_model_file = std::stoi(num_vars_string); + nvars_model_file = std::stoi(num_vars_string); // now read in the names of the variables - std::vector varnames_stored; for (int n = 0; n < nvars_model_file; n++) { getline(mf, line); std::string var_string = line.substr(line.find("#")+1, line.length()); varnames_stored.push_back(model_string::ltrim(model_string::rtrim(var_string))); } } else { - // new format + // new format. example: + """ + # radius density temperature + 0.1 3 4.0 + 0.2 4 3.9 + 0.3 4 3.7 + """ + + // remove the leading "#" + if (!line.empty() && line[0] == '#') { + line = line.substr(1); + } - + std::istringstream iss(line); + + // in this format, the first column name is the coordinate -- + // we can skip that + + std::string word; + iss >> word; + + while (iss >> word) { + varnames_stored.push_back(word); + } + + nvars_model_file = static_cast(varnames_stored.size()); } From a2b301af43a432c108af3ec7f28eb8fb2066e3da Mon Sep 17 00:00:00 2001 From: Melissa Rasmussen Date: Thu, 5 Feb 2026 21:17:20 -0500 Subject: [PATCH 03/10] Fix data input to generalize between old and new formats --- read_model.H | 31 ++++++++++++++++++------------- 1 file changed, 18 insertions(+), 13 deletions(-) diff --git a/read_model.H b/read_model.H index 9eb37ea..11457c2 100644 --- a/read_model.H +++ b/read_model.H @@ -162,8 +162,6 @@ read_file(const std::string filename, initial_model_t& initial_model) { getline(mf, line); - // variables: - // initial_model.npts int nvars_model_file; std::vector varnames_stored; @@ -180,12 +178,8 @@ read_file(const std::string filename, initial_model_t& initial_model) { 0.3 4 3.7 """ - std::string npts_string = line.substr(line.find("=")+1, line.length()); - initial_model.npts = std::stoi(npts_string); - - if (initial_model.npts > NPTS_MODEL) { - amrex::Error("Error: model has more than NPTS_MODEL points,"); - } + // ignore the first line, and detect the number of points + // automatically below // next line tells use the number of variables @@ -240,11 +234,17 @@ read_file(const std::string filename, initial_model_t& initial_model) { amrex::Vector vars_stored; vars_stored.resize(nvars_model_file); - for (int i = 0; i < initial_model.npts; i++) { - mf >> initial_model.r(i); + int i{0}; + while (getline(mf, line)) { + if (i > NPTS_MODEL) { + amrex::Error("Error: model has more than NPTS_MODEL points. Increase MAX_NPTS_MODEL"); + } + + std::istringstream iss(line); + iss >> initial_model.r(i); for (int j = 0; j < nvars_model_file; j++) { - mf >> vars_stored[j]; + iss >> vars_stored[j]; } for (int j = 0; j < model::nvar; j++) { @@ -290,9 +290,14 @@ read_file(const std::string filename, initial_model_t& initial_model) { amrex::Print() << Font::Bold << FGColor::Yellow << "[WARNING] variable not found: " << varnames_stored[j] << ResetDisplay << std::endl; } - } // end loop over nvars_model_file + } // end loop over variables in line - } + i++; + } // end of loop over lines in model file + + initial_model.npts = i; + + mf.close(); } #endif From 37901d5bc0c7b2f598f84ca6286c179cf9539fdc Mon Sep 17 00:00:00 2001 From: Melissa Rasmussen Date: Thu, 5 Feb 2026 21:25:06 -0500 Subject: [PATCH 04/10] Debug errors; backward compatibility checked --- read_model.H | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/read_model.H b/read_model.H index 11457c2..c234794 100644 --- a/read_model.H +++ b/read_model.H @@ -167,7 +167,7 @@ read_file(const std::string filename, initial_model_t& initial_model) { if (line.find("npts") != std::string::npos) { // old format. example: - """ + /* # npts = 3 # number of variables = 2 # radius @@ -176,7 +176,7 @@ read_file(const std::string filename, initial_model_t& initial_model) { 0.1 3 4.0 0.2 4 3.9 0.3 4 3.7 - """ + */ // ignore the first line, and detect the number of points // automatically below @@ -185,7 +185,7 @@ read_file(const std::string filename, initial_model_t& initial_model) { getline(mf, line); - if (line.find("num of variables") == std::string::npos) { + if (line.find("num of variables") == std::string::npos && line.find("number of variables") == std::string::npos) { amrex::Error("Error: missing the number of variables header"); } @@ -201,12 +201,12 @@ read_file(const std::string filename, initial_model_t& initial_model) { } } else { // new format. example: - """ + /* # radius density temperature 0.1 3 4.0 0.2 4 3.9 0.3 4 3.7 - """ + */ // remove the leading "#" if (!line.empty() && line[0] == '#') { From cfc56b0831080ece52289d0a50c64665b5d7f5dd Mon Sep 17 00:00:00 2001 From: Melissa Rasmussen Date: Thu, 12 Feb 2026 14:38:16 -0500 Subject: [PATCH 05/10] make identification of old/new file format more broad. Fixes the lagrangian_planar case --- read_model.H | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/read_model.H b/read_model.H index c234794..7a26d3b 100644 --- a/read_model.H +++ b/read_model.H @@ -165,7 +165,7 @@ read_file(const std::string filename, initial_model_t& initial_model) { int nvars_model_file; std::vector varnames_stored; - if (line.find("npts") != std::string::npos) { + if (line.find("npts") != std::string::npos || line.find("number of points") != std::string::npos) { // old format. example: /* # npts = 3 From 8334d9fbc458db2b3615b5dd15cb951f14653940 Mon Sep 17 00:00:00 2001 From: Melissa Rasmussen Date: Thu, 12 Feb 2026 15:10:38 -0500 Subject: [PATCH 06/10] bring model file (and the script that generates it) into compliance with Castro old format --- massive_star/15m_500_sec.aprox19.dat | 55 ++++++++++++++-------------- massive_star/convert_21_to_19.py | 5 +-- 2 files changed, 28 insertions(+), 32 deletions(-) diff --git a/massive_star/15m_500_sec.aprox19.dat b/massive_star/15m_500_sec.aprox19.dat index e92e4f6..47a5149 100644 --- a/massive_star/15m_500_sec.aprox19.dat +++ b/massive_star/15m_500_sec.aprox19.dat @@ -1,32 +1,32 @@ # npts = 3611 # number of variables = 27 -density -temperature -ye -velx -pressure -entr -abar -enuc -n -H1 -p -He3 -He4 -C12 -N14 -O16 -Ne20 -Mg24 -Si28 -S32 -Ar36 -Ca40 -Ti44 -Cr48 -Fe52 -Fe54 -Ni56 +# density +# temperature +# ye +# velx +# pressure +# entr +# abar +# enuc +# n +# H1 +# p +# He3 +# He4 +# C12 +# N14 +# O16 +# Ne20 +# Mg24 +# Si28 +# S32 +# Ar36 +# Ca40 +# Ti44 +# Cr48 +# Fe52 +# Fe54 +# Ni56 337537.71121300000232 1808953118.9500000477 6584485414.5799999237 0.46354931956700001772 -119.34958673799999929 1.0507722924199999529e+27 0.85488876350700004902 51.959996281300000476 10717328463399999488 0.0012822561731099999232 0 8.734243857490000473e-08 0 0.00055562984419200002559 3.5976537168900001507e-09 0 7.9770873314700001431e-10 2.5784726164900001073e-12 2.784609616930000066e-11 1.0159525520600000283e-09 1.5988627581300000157e-10 1.8389748363700000301e-11 3.5904618100499999432e-12 3.0906166241799997709e-14 2.7972765513299999597e-14 2.4717986882499998539e-14 2.603937801570000207e-06 0.99815941707886168821 594047.20889500004705 1808775668.3699998856 6584256849.8100004196 0.46354956795199997455 -150.37372296300000585 1.0506333147799999863e+27 0.85488120579199999849 51.961420944099998565 -3028276060670000128 0.001281773310189999928 0 8.7287310247000001032e-08 0 0.00055539992557500003894 3.594294609229999968e-09 0 7.969471086430000436e-10 2.5755498214499998767e-12 2.7817658838999998508e-11 1.0150651159099999464e-09 1.5974328974399999034e-10 1.8372255479600000789e-11 3.5869898665700002016e-12 3.0871376680500002269e-14 2.7942324735400000559e-14 2.4692991233799999298e-14 2.6028843046100000086e-06 0.99816013097389177933 832856.79917100002058 1808546203.4000000954 6584382058.4099998474 0.46354961633300001411 -189.46412342799999351 1.0504678768300000408e+27 0.85492014537299998445 51.961226871800000993 27715722773700001792 0.0012817750373099998963 0 8.7383507065000006563e-08 0 0.00055570183621600000631 3.598040236700000117e-09 0 7.9785030440799997164e-10 2.5789263096800000188e-12 2.7854509254399999653e-11 1.0164074134999999535e-09 1.5997014328799999378e-10 1.8400719411899998846e-11 3.5929843078699998736e-12 3.0929638242799998024e-14 2.7998417638099998454e-14 2.4744858893900000256e-14 2.6054742517499999882e-06 0.99815982464369179539 @@ -3638,4 +3638,3 @@ Ni56 63436906460500 2.1977853600099999171e-09 3168.5209801699998025 0.82242463736199999946 2015.2893995800000084 419.88098907900001677 20.322542062800000195 1.3700919410100000917 2.7360769874299997925e-28 1.0035163404800000013e-99 0.64495034955999996917 1.0035163404800000013e-99 1.9571546147100000226e-05 0.33498340974800000502 0.0020667958887200001128 0.0040157717779699996299 0.0079674516039400008344 0.0021288213024800001737 0.00081081192470599996319 0.00089476020801399999173 0.00044591876249999999778 8.8048009653300004476e-05 7.8728240015199998383e-05 1.4923075529999999622e-68 1.0035163404800000013e-99 1.0035163404800000013e-99 8.9681194949700002985e-05 0.0014598802325300000177 63436906470600 2.1977853541799999406e-09 3168.5209792800001196 0.82242463736399995522 2015.2768581400000585 419.88098784499999283 20.322542064400000328 1.3700919410100000917 2.7360769802199999217e-28 1.0035163404800000013e-99 0.64495034956499996959 1.0035163404800000013e-99 1.9571546147700001096e-05 0.33498340974300000461 0.002066795888779999861 0.0040157717778000001987 0.0079674516040600003308 0.0021288213024800001737 0.00081081192470599996319 0.00089476020801399999173 0.00044591876249999999778 8.8048009653300004476e-05 7.8728240015199998383e-05 1.4923075527700000884e-68 1.0035163404800000013e-99 1.0035163404800000013e-99 8.9681194949700002985e-05 0.0014598802325300000177 63436906475700 2.1977853522400000151e-09 3168.520978980000109 0.822424637366000022 2015.2768581400000585 419.88098743400001922 20.322542065000000377 1.3700919410000000909 2.7360769778199998534e-28 1.0035163404800000013e-99 0.64495034956899999212 1.0035163404800000013e-99 1.9571546148100000547e-05 0.33498340974000001546 0.0020667958888199999823 0.0040157717776700000215 0.0079674516041400005734 0.0021288213024800001737 0.00081081192470599996319 0.00089476020801399999173 0.00044591876249999999778 8.8048009653300004476e-05 7.8728240015199998383e-05 1.4923075525900000887e-68 1.0035163404800000013e-99 1.0035163404800000013e-99 8.9681194949700002985e-05 0.0014598802325300000177 -# conversion of 15m_500_sec.txt to aprox19 nuclei via convert_21_to_19.py diff --git a/massive_star/convert_21_to_19.py b/massive_star/convert_21_to_19.py index 2001018..c8e4369 100644 --- a/massive_star/convert_21_to_19.py +++ b/massive_star/convert_21_to_19.py @@ -95,11 +95,8 @@ for name in names_new: if name == "r": continue - f.write(f"{name}\n") + f.write(f"# {name}\n") for irow in range(data_new.shape[0]): l = [f"{q:30.20g}" for q in data_new[irow, :]] f.write(" ".join(l) + "\n") - f.write(f"# conversion of {file} to aprox19 nuclei via convert_21_to_19.py\n") - - From 554ea3cb57a7819f63e52aefaa8d38a16a2121d1 Mon Sep 17 00:00:00 2001 From: Melissa Rasmussen Date: Fri, 13 Feb 2026 13:21:31 -0500 Subject: [PATCH 07/10] fix mass calculation to work for both ascending and descending radius values --- massive_star/init_1d.H | 16 ++++++++++------ 1 file changed, 10 insertions(+), 6 deletions(-) diff --git a/massive_star/init_1d.H b/massive_star/init_1d.H index a84bf6b..29612bd 100644 --- a/massive_star/init_1d.H +++ b/massive_star/init_1d.H @@ -101,8 +101,6 @@ AMREX_INLINE void init_1d() { // are mapping onto, and then we want to force it into HSE on that // mesh. - std::cout << "here " << problem_rp::nx << " " << NPTS_MODEL << std::endl; - if (problem_rp::nx > NPTS_MODEL) { amrex::Error("too many zones requested -- increase NPTS_MODEL"); } @@ -144,23 +142,29 @@ AMREX_INLINE void init_1d() { of.close(); + // determine whether the radial coordinates are ascending or descending. + bool ascending = initial_model.r(initial_model.npts-1) > initial_model.r(0); + // compute the mass of the initial model -- the data is // non-uniform, and appears to be node-centered finite difference. // We'll just do a very simple trapezoid rule. // the first node is not at r = 0, so just use a rectangle rule to get us started - Real M_total = (4.0_rt / 3.0_rt) * M_PI * std::pow(initial_model.r(0), 3) * - initial_model.state(0, model::idens); + Real smallest_r_index = ascending ? 0 : initial_model.npts-1; + Real M_total = (4.0_rt / 3.0_rt) * M_PI * std::pow(initial_model.r(smallest_r_index), 3) * + initial_model.state(smallest_r_index, model::idens); - for (int i = 1; i < initial_model.npts-1; ++i) { + for (int i = 1; i < initial_model.npts; ++i) { // integrate 4pi r**2 rho using trapezoid Real rl = initial_model.r(i-1); Real rr = initial_model.r(i); - M_total += 4.0_rt * M_PI * 0.5_rt * (rr - rl) * + Real diff = ascending ? rr - rl : rl - rr; + + M_total += 4.0_rt * M_PI * 0.5_rt * diff * (rl * rl * initial_model.state(i-1, model::idens) + rr * rr * initial_model.state(i, model::idens)); } From eb044b08c72fb348dea6e4d574618a6a4c2b0b45 Mon Sep 17 00:00:00 2001 From: Melissa Rasmussen Date: Fri, 13 Feb 2026 14:45:11 -0500 Subject: [PATCH 08/10] Move ascending/descending fix into read_model --- massive_star/init_1d.H | 14 +++++--------- read_model.H | 19 +++++++++++++++++++ 2 files changed, 24 insertions(+), 9 deletions(-) diff --git a/massive_star/init_1d.H b/massive_star/init_1d.H index 29612bd..96c5d6d 100644 --- a/massive_star/init_1d.H +++ b/massive_star/init_1d.H @@ -129,6 +129,7 @@ AMREX_INLINE void init_1d() { std::ofstream of; of.open("model.orig"); + // note the rows will be in reversed order if input file is descending in radius of << "# initial model as read in" << std::endl; @@ -142,18 +143,14 @@ AMREX_INLINE void init_1d() { of.close(); - // determine whether the radial coordinates are ascending or descending. - bool ascending = initial_model.r(initial_model.npts-1) > initial_model.r(0); - // compute the mass of the initial model -- the data is // non-uniform, and appears to be node-centered finite difference. // We'll just do a very simple trapezoid rule. // the first node is not at r = 0, so just use a rectangle rule to get us started - Real smallest_r_index = ascending ? 0 : initial_model.npts-1; - Real M_total = (4.0_rt / 3.0_rt) * M_PI * std::pow(initial_model.r(smallest_r_index), 3) * - initial_model.state(smallest_r_index, model::idens); + Real M_total = (4.0_rt / 3.0_rt) * M_PI * std::pow(initial_model.r(0), 3) * + initial_model.state(0, model::idens); for (int i = 1; i < initial_model.npts; ++i) { @@ -162,13 +159,12 @@ AMREX_INLINE void init_1d() { Real rl = initial_model.r(i-1); Real rr = initial_model.r(i); - Real diff = ascending ? rr - rl : rl - rr; - - M_total += 4.0_rt * M_PI * 0.5_rt * diff * + M_total += 4.0_rt * M_PI * 0.5_rt * (rr - rl) * (rl * rl * initial_model.state(i-1, model::idens) + rr * rr * initial_model.state(i, model::idens)); } + Real smallest_r_index = ascending ? 0 : initial_model.npts-1; std::cout << "total mass as read from original model = " << M_total << std::endl; // put the model onto our new uniform grid diff --git a/read_model.H b/read_model.H index 7a26d3b..6ceffbd 100644 --- a/read_model.H +++ b/read_model.H @@ -298,6 +298,25 @@ read_file(const std::string filename, initial_model_t& initial_model) { initial_model.npts = i; mf.close(); + + // ensure model is in order of ascending radius, not descending + if (initial_model.r(initial_model.npts-1) > initial_model.r(0)) { + amrex::Array2D state_copy; + amrex::Array2D r_copy; + + for (int i = 0; i < initial_model.npts; ++i){ // copy + r_copy(i) = initial_model.r(i); + for (int j = 0; j < model::nvar; ++j) { + state_copy(i,j) = initial_model.state(i,j); + } + } + for (int i = 0; i < initial_model.npts; ++i){ // reverse + initial_model.r(i) = r_copy(initial_model.npts-1-i); + for (int j = 0; j < model::nvar; ++j) { + initial_model.state(i,j) = state_copy(initial_model.npts-1-i,j); + } + } + } } #endif From 9d04abc595b2b086b6b08531afd5b131faa3449e Mon Sep 17 00:00:00 2001 From: Melissa Rasmussen Date: Fri, 13 Feb 2026 15:15:05 -0500 Subject: [PATCH 09/10] debugging: reverse init model causes segfault --- massive_star/init_1d.H | 1 - read_model.H | 9 ++++++--- 2 files changed, 6 insertions(+), 4 deletions(-) diff --git a/massive_star/init_1d.H b/massive_star/init_1d.H index 96c5d6d..299a0ac 100644 --- a/massive_star/init_1d.H +++ b/massive_star/init_1d.H @@ -164,7 +164,6 @@ AMREX_INLINE void init_1d() { rr * rr * initial_model.state(i, model::idens)); } - Real smallest_r_index = ascending ? 0 : initial_model.npts-1; std::cout << "total mass as read from original model = " << M_total << std::endl; // put the model onto our new uniform grid diff --git a/read_model.H b/read_model.H index 6ceffbd..cdf9855 100644 --- a/read_model.H +++ b/read_model.H @@ -300,20 +300,23 @@ read_file(const std::string filename, initial_model_t& initial_model) { mf.close(); // ensure model is in order of ascending radius, not descending - if (initial_model.r(initial_model.npts-1) > initial_model.r(0)) { + if (initial_model.r(initial_model.npts-1) < initial_model.r(0)) { amrex::Array2D state_copy; - amrex::Array2D r_copy; + amrex::Array1D r_copy; for (int i = 0; i < initial_model.npts; ++i){ // copy r_copy(i) = initial_model.r(i); for (int j = 0; j < model::nvar; ++j) { state_copy(i,j) = initial_model.state(i,j); + std::cout << state_copy(i,j) << " "; } } + std::cout << std::endl << "compared to:" << std::endl; for (int i = 0; i < initial_model.npts; ++i){ // reverse initial_model.r(i) = r_copy(initial_model.npts-1-i); for (int j = 0; j < model::nvar; ++j) { - initial_model.state(i,j) = state_copy(initial_model.npts-1-i,j); + //std::cout << state_copy(i,j) << " "; //THIS LINE CAUSES SEGFAULT IF UNCOMMENTED + //initial_model.state(i,j) = state_copy(initial_model.npts-1-i,j); } } } From 012537febeeaed880d033fb9caa014e4ed73a0a5 Mon Sep 17 00:00:00 2001 From: Melissa Rasmussen Date: Wed, 18 Feb 2026 16:21:25 -0500 Subject: [PATCH 10/10] Make swap memory-efficient, avoiding segfault --- read_model.H | 50 ++++++++++++++++++++++++-------------------------- 1 file changed, 24 insertions(+), 26 deletions(-) diff --git a/read_model.H b/read_model.H index cdf9855..d7c5f0b 100644 --- a/read_model.H +++ b/read_model.H @@ -234,21 +234,21 @@ read_file(const std::string filename, initial_model_t& initial_model) { amrex::Vector vars_stored; vars_stored.resize(nvars_model_file); - int i{0}; + int line_num{0}; while (getline(mf, line)) { - if (i > NPTS_MODEL) { + if (line_num > NPTS_MODEL) { amrex::Error("Error: model has more than NPTS_MODEL points. Increase MAX_NPTS_MODEL"); } std::istringstream iss(line); - iss >> initial_model.r(i); + iss >> initial_model.r(line_num); for (int j = 0; j < nvars_model_file; j++) { iss >> vars_stored[j]; } for (int j = 0; j < model::nvar; j++) { - initial_model.state(i,j) = 0.0_rt; + initial_model.state(line_num,j) = 0.0_rt; } @@ -259,25 +259,25 @@ read_file(const std::string filename, initial_model_t& initial_model) { bool found_model = false; if (varnames_stored[j] == "density") { - initial_model.state(i,model::idens) = vars_stored[j]; + initial_model.state(line_num,model::idens) = vars_stored[j]; found_model = true; } else if (varnames_stored[j] == "temperature") { - initial_model.state(i,model::itemp) = vars_stored[j]; + initial_model.state(line_num,model::itemp) = vars_stored[j]; found_model = true; } else if (varnames_stored[j] == "pressure") { - initial_model.state(i,model::ipres) = vars_stored[j]; + initial_model.state(line_num,model::ipres) = vars_stored[j]; found_model = true; } else if (varnames_stored[j] == "ye" || varnames_stored[j] == "Ye") { - initial_model.state(i,model::iyef) = vars_stored[j]; + initial_model.state(line_num,model::iyef) = vars_stored[j]; found_model = true; } else { for (int comp = 0; comp < NumSpec; comp++) { if (varnames_stored[j] == spec_names_cxx[comp] || varnames_stored[j] == short_spec_names_cxx[comp]) { - initial_model.state(i,model::ispec+comp) = vars_stored[j]; + initial_model.state(line_num,model::ispec+comp) = vars_stored[j]; found_model = true; break; } @@ -286,37 +286,35 @@ read_file(const std::string filename, initial_model_t& initial_model) { // yell if we didn't find the current variable - if (!found_model && i == 0) { + if (!found_model && line_num == 0) { amrex::Print() << Font::Bold << FGColor::Yellow << "[WARNING] variable not found: " << varnames_stored[j] << ResetDisplay << std::endl; } } // end loop over variables in line - i++; + line_num++; } // end of loop over lines in model file - initial_model.npts = i; + initial_model.npts = line_num; mf.close(); // ensure model is in order of ascending radius, not descending if (initial_model.r(initial_model.npts-1) < initial_model.r(0)) { - amrex::Array2D state_copy; - amrex::Array1D r_copy; + std::cout << "in if" << std::endl; + amrex::Real swap_var; + + for (int i = 0; i < int(std::floor(initial_model.npts/2)); ++i){ + int swap_row = initial_model.npts-1-i; + + swap_var = initial_model.r(i); + initial_model.r(i) = initial_model.r(swap_row); + initial_model.r(swap_row) = swap_var; - for (int i = 0; i < initial_model.npts; ++i){ // copy - r_copy(i) = initial_model.r(i); - for (int j = 0; j < model::nvar; ++j) { - state_copy(i,j) = initial_model.state(i,j); - std::cout << state_copy(i,j) << " "; - } - } - std::cout << std::endl << "compared to:" << std::endl; - for (int i = 0; i < initial_model.npts; ++i){ // reverse - initial_model.r(i) = r_copy(initial_model.npts-1-i); for (int j = 0; j < model::nvar; ++j) { - //std::cout << state_copy(i,j) << " "; //THIS LINE CAUSES SEGFAULT IF UNCOMMENTED - //initial_model.state(i,j) = state_copy(initial_model.npts-1-i,j); + swap_var = initial_model.state(i,j); + initial_model.state(i,j) = initial_model.state(swap_row,j); + initial_model.state(swap_row,j) = swap_var; } } }