mirror of
https://github.com/cosmo-sims/monofonIC.git
synced 2024-09-19 17:03:45 +02:00
added Bodo Schwabe's work in progress NyX plugin implementation
This commit is contained in:
parent
326bfb7b9e
commit
210bde40fb
3 changed files with 610 additions and 0 deletions
|
@ -273,3 +273,11 @@ if(ENABLE_GENERICIO)
|
||||||
target_link_libraries(${PRGNAME} PRIVATE genericio::genericio_mpi)
|
target_link_libraries(${PRGNAME} PRIVATE genericio::genericio_mpi)
|
||||||
target_compile_definitions(${PRGNAME} PRIVATE "ENABLE_GENERICIO")
|
target_compile_definitions(${PRGNAME} PRIVATE "ENABLE_GENERICIO")
|
||||||
endif()
|
endif()
|
||||||
|
|
||||||
|
option(ENABLE_AMReX "Enable AMReX support (as a submodule)." OFF)
|
||||||
|
|
||||||
|
if(ENABLE_AMReX)
|
||||||
|
include(${CMAKE_CURRENT_SOURCE_DIR}/external/AMReX.cmake)
|
||||||
|
target_link_libraries( ${PRGNAME} PRIVATE AMReX::amrex AMReX::Flags_CXX )
|
||||||
|
target_compile_definitions(${PRGNAME} PRIVATE "ENABLE_AMREX")
|
||||||
|
endif()
|
17
external/AMReX.cmake
vendored
Normal file
17
external/AMReX.cmake
vendored
Normal file
|
@ -0,0 +1,17 @@
|
||||||
|
cmake_minimum_required(VERSION 3.11)
|
||||||
|
include(FetchContent)
|
||||||
|
FetchContent_Declare(
|
||||||
|
amrex
|
||||||
|
GIT_REPOSITORY https://github.com/AMReX-Codes/amrex.git
|
||||||
|
GIT_TAG development
|
||||||
|
GIT_SHALLOW YES
|
||||||
|
GIT_PROGRESS TRUE
|
||||||
|
USES_TERMINAL_DOWNLOAD TRUE # <---- this is needed only for Ninja
|
||||||
|
)
|
||||||
|
|
||||||
|
FetchContent_GetProperties(amrex)
|
||||||
|
if(NOT amrex)
|
||||||
|
set(FETCHCONTENT_QUIET OFF)
|
||||||
|
FetchContent_Populate(amrex)
|
||||||
|
add_subdirectory(${amrex_SOURCE_DIR} ${amrex_BINARY_DIR})
|
||||||
|
endif()
|
585
src/plugins/output_nyx.cc
Normal file
585
src/plugins/output_nyx.cc
Normal file
|
@ -0,0 +1,585 @@
|
||||||
|
// This file is part of monofonIC (MUSIC2)
|
||||||
|
// A software package to generate ICs for cosmological simulations
|
||||||
|
// Copyright (C) 2020 by Oliver Hahn
|
||||||
|
// Copyright (C) 2012 by Jan Frederik Engels
|
||||||
|
// Copyright (C) 2020 by Bodo Schwabe
|
||||||
|
//
|
||||||
|
// monofonIC is free software: you can redistribute it and/or modify
|
||||||
|
// it under the terms of the GNU General Public License as published by
|
||||||
|
// the Free Software Foundation, either version 3 of the License, or
|
||||||
|
// (at your option) any later version.
|
||||||
|
//
|
||||||
|
// monofonIC is distributed in the hope that it will be useful,
|
||||||
|
// but WITHOUT ANY WARRANTY; without even the implied warranty of
|
||||||
|
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
||||||
|
// GNU General Public License for more details.
|
||||||
|
//
|
||||||
|
// You should have received a copy of the GNU General Public License
|
||||||
|
// along with this program. If not, see <http://www.gnu.org/licenses/>.
|
||||||
|
/*
|
||||||
|
|
||||||
|
output_nyx.cc - This file is part of MUSIC -
|
||||||
|
a code to generate multi-scale initial conditions
|
||||||
|
for cosmological simulations
|
||||||
|
|
||||||
|
Copyright (C) 2010 Oliver Hahn
|
||||||
|
Copyright (C) 2012 Jan Frederik Engels
|
||||||
|
|
||||||
|
*/
|
||||||
|
|
||||||
|
#ifdef ENABLE_AMREX
|
||||||
|
|
||||||
|
#include <output_plugin.hh>
|
||||||
|
|
||||||
|
#include <AMReX_VisMF.H>
|
||||||
|
#include <AMReX_Box.H>
|
||||||
|
#include <AMReX_RealBox.H>
|
||||||
|
#include <AMReX_ParallelDescriptor.H>
|
||||||
|
#include <AMReX_Utility.H>
|
||||||
|
#include <AMReX_Array.H>
|
||||||
|
#include <AMReX_BoxArray.H>
|
||||||
|
#include <AMReX_FabArray.H>
|
||||||
|
#include <AMReX_MultiFab.H>
|
||||||
|
|
||||||
|
#define MAX_GRID_SIZE 32
|
||||||
|
#define BL_SPACEDIM 3
|
||||||
|
|
||||||
|
class nyx_output_plugin : public output_plugin
|
||||||
|
{
|
||||||
|
protected:
|
||||||
|
|
||||||
|
struct patch_header{
|
||||||
|
int component_rank;
|
||||||
|
size_t component_size;
|
||||||
|
std::vector<int> dimensions;
|
||||||
|
int rank;
|
||||||
|
std::vector<int> top_grid_dims;
|
||||||
|
std::vector<int> top_grid_end;
|
||||||
|
std::vector<int> top_grid_start;
|
||||||
|
};
|
||||||
|
|
||||||
|
struct sim_header{
|
||||||
|
std::vector<int> dimensions;
|
||||||
|
std::vector<int> offset;
|
||||||
|
float a_start;
|
||||||
|
float dx;
|
||||||
|
float h0;
|
||||||
|
float omega_b;
|
||||||
|
float omega_m;
|
||||||
|
float omega_v;
|
||||||
|
float vfact;
|
||||||
|
float boxlength;
|
||||||
|
int particle_idx;
|
||||||
|
};
|
||||||
|
|
||||||
|
int n_data_items;
|
||||||
|
std::vector<std::string> field_name;
|
||||||
|
int f_lev;
|
||||||
|
int gridp;
|
||||||
|
uint32_t levelmin_;
|
||||||
|
uint32_t levelmax_;
|
||||||
|
real_t lunit_, vunit_, munit_;
|
||||||
|
|
||||||
|
std::vector<amrex::MultiFab*> mfs;
|
||||||
|
|
||||||
|
std::vector<amrex::BoxArray> boxarrays;
|
||||||
|
std::vector<amrex::Box> boxes;
|
||||||
|
|
||||||
|
sim_header the_sim_header;
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
// void dump_grid_data(int comp, std::string fieldname, const grid_hierarchy& gh, double factor = 1.0, double add = 0.0 )
|
||||||
|
// {
|
||||||
|
// std::cout << fieldname << " is dumped... to mf index " << comp << std::endl;
|
||||||
|
|
||||||
|
// //FIXME adapt for multiple levels!
|
||||||
|
// for(int mlevel=levelmin_; mlevel<=levelmax_; ++mlevel )
|
||||||
|
// {
|
||||||
|
// int blevel = mlevel-levelmin_;
|
||||||
|
|
||||||
|
// std::vector<int> ng;
|
||||||
|
// ng.push_back( gh.get_grid(mlevel)->size(0) );
|
||||||
|
// ng.push_back( gh.get_grid(mlevel)->size(1) );
|
||||||
|
// ng.push_back( gh.get_grid(mlevel)->size(2) );
|
||||||
|
|
||||||
|
// std::cout << ng[0] << " " << ng[1] << " " << ng[2] << std::endl;
|
||||||
|
|
||||||
|
// //write data to mf
|
||||||
|
// for(MFIter mfi(*(mfs[blevel])); mfi.isValid(); ++mfi) {
|
||||||
|
// FArrayBox &myFab = (*(mfs[blevel]))[mfi];
|
||||||
|
// const Box& box = mfi.validbox();
|
||||||
|
// const int *fab_lo = box.loVect();
|
||||||
|
// const int *fab_hi = box.hiVect();
|
||||||
|
|
||||||
|
// int mk = fab_lo[2] - boxes[blevel].smallEnd()[2];
|
||||||
|
// #ifdef OMP
|
||||||
|
// #pragma omp parallel for default(shared)
|
||||||
|
// #endif
|
||||||
|
// for (int k = fab_lo[2]; k <= fab_hi[2]; k++, mk++) {
|
||||||
|
|
||||||
|
// int mj = fab_lo[1] - boxes[blevel].smallEnd()[1];
|
||||||
|
// for (int j = fab_lo[1]; j <= fab_hi[1]; j++, mj++) {
|
||||||
|
|
||||||
|
// int mi = fab_lo[0] - boxes[blevel].smallEnd()[0];
|
||||||
|
// for (int i = fab_lo[0]; i <= fab_hi[0]; i++, mi++) {
|
||||||
|
// if (mi>=ng[0])
|
||||||
|
// std::cout << "mi (" << mi << ") too large " << ng[0] << std::endl;
|
||||||
|
// if (mj>=ng[1])
|
||||||
|
// std::cout << "mj (" << mj << ") too large " << ng[1] << std::endl;
|
||||||
|
// if (mk>=ng[2])
|
||||||
|
// std::cout << "mk (" << mk << ") too large " << ng[2] << std::endl;
|
||||||
|
|
||||||
|
// IntVect iv(i,j,k);
|
||||||
|
// double data = ( add + (*gh.get_grid(mlevel))(mi,mj,mk) )*factor;
|
||||||
|
// int idx = myFab.box().index(iv);
|
||||||
|
// myFab.dataPtr(comp)[idx] = data;
|
||||||
|
|
||||||
|
// }
|
||||||
|
// }
|
||||||
|
// }
|
||||||
|
// }
|
||||||
|
// }
|
||||||
|
// }
|
||||||
|
|
||||||
|
public:
|
||||||
|
//constructor
|
||||||
|
explicit nyx_output_plugin( config_file& cf, std::unique_ptr<cosmology::calculator> &pcc )
|
||||||
|
: output_plugin(cf, pcc, "NYX/AMREX")
|
||||||
|
{
|
||||||
|
int argc=0;
|
||||||
|
char **argv;
|
||||||
|
amrex::Initialize(argc,argv);
|
||||||
|
|
||||||
|
bool bhave_hydro = cf_.get_value_safe<bool>("setup", "baryons", false);
|
||||||
|
|
||||||
|
if (bhave_hydro)
|
||||||
|
n_data_items = 10;
|
||||||
|
else
|
||||||
|
n_data_items = 6;
|
||||||
|
|
||||||
|
field_name.resize(n_data_items);
|
||||||
|
if (bhave_hydro)
|
||||||
|
{
|
||||||
|
field_name[0] = "baryon_density";
|
||||||
|
field_name[1] = "baryon_vel_x";
|
||||||
|
field_name[2] = "baryon_vel_y";
|
||||||
|
field_name[3] = "baryon_vel_z";
|
||||||
|
field_name[4] = "dm_pos_x";
|
||||||
|
field_name[5] = "dm_pos_y";
|
||||||
|
field_name[6] = "dm_pos_z";
|
||||||
|
field_name[7] = "dm_vel_x";
|
||||||
|
field_name[8] = "dm_vel_y";
|
||||||
|
field_name[9] = "dm_vel_z";
|
||||||
|
the_sim_header.particle_idx = 4;
|
||||||
|
}
|
||||||
|
else
|
||||||
|
{
|
||||||
|
field_name[0] = "dm_pos_x";
|
||||||
|
field_name[1] = "dm_pos_y";
|
||||||
|
field_name[2] = "dm_pos_z";
|
||||||
|
field_name[3] = "dm_vel_x";
|
||||||
|
field_name[4] = "dm_vel_y";
|
||||||
|
field_name[5] = "dm_vel_z";
|
||||||
|
the_sim_header.particle_idx = 0;
|
||||||
|
}
|
||||||
|
|
||||||
|
uint32_t ngrid = cf_.get_value<int>("setup", "GridRes");
|
||||||
|
levelmin_ = uint32_t(std::log2(double(ngrid)) + 1e-6);
|
||||||
|
|
||||||
|
//for future multilevel simulations
|
||||||
|
levelmax_ = uint32_t(std::log2(double(ngrid)) + 1e-6);
|
||||||
|
double offx_[]={0.0};
|
||||||
|
double offy_[]={0.0};
|
||||||
|
double offz_[]={0.0};
|
||||||
|
int sizex_[] = {int(ngrid)};
|
||||||
|
int sizey_[] = {int(ngrid)};
|
||||||
|
int sizez_[] = {int(ngrid)};
|
||||||
|
|
||||||
|
f_lev = levelmax_-levelmin_;
|
||||||
|
std::cout << f_lev+1 << " level" << std::endl;
|
||||||
|
|
||||||
|
mfs.resize(f_lev+1);
|
||||||
|
|
||||||
|
amrex::Vector<int> pmap(2);
|
||||||
|
pmap[0]=0;
|
||||||
|
pmap[1]=0;
|
||||||
|
|
||||||
|
gridp = 1<<levelmin_; //probably ngrid
|
||||||
|
|
||||||
|
double off[] = {0, 0, 0};
|
||||||
|
|
||||||
|
//at first we do this only for the topgrid...
|
||||||
|
for(int lev = 0; lev <= f_lev; lev++)
|
||||||
|
{
|
||||||
|
amrex::BoxArray domainBoxArray(1);
|
||||||
|
|
||||||
|
off[0] += offx_[lev];
|
||||||
|
off[1] += offy_[lev];
|
||||||
|
off[2] += offz_[lev];
|
||||||
|
|
||||||
|
for (int asdf = 0; asdf < 3; asdf++)
|
||||||
|
off[asdf] *= 2;
|
||||||
|
|
||||||
|
amrex::IntVect pdLo(off[0],
|
||||||
|
off[1],
|
||||||
|
off[2]);
|
||||||
|
amrex::IntVect pdHi(off[0]+sizex_[lev]-1,
|
||||||
|
off[1]+sizey_[lev]-1,
|
||||||
|
off[2]+sizez_[lev]-1);
|
||||||
|
|
||||||
|
std::cout << pdLo << std::endl;
|
||||||
|
std::cout << pdHi << std::endl;
|
||||||
|
|
||||||
|
// Start with a probDomain
|
||||||
|
amrex::Box probDomain(pdLo,pdHi);
|
||||||
|
|
||||||
|
// We just have one box since we don't use mpi.
|
||||||
|
domainBoxArray.set(0, probDomain);
|
||||||
|
domainBoxArray.maxSize(32);
|
||||||
|
pmap.resize(domainBoxArray.size(),0);
|
||||||
|
|
||||||
|
amrex::DistributionMapping domainDistMap(pmap);
|
||||||
|
|
||||||
|
boxarrays.push_back(domainBoxArray);
|
||||||
|
boxes.push_back(probDomain);
|
||||||
|
|
||||||
|
int ngrow(0);
|
||||||
|
mfs[lev] = new amrex::MultiFab(domainBoxArray, domainDistMap, n_data_items, ngrow);
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
bool haveblockingfactor = cf.contains_key( "setup", "blocking_factor");
|
||||||
|
|
||||||
|
if( !haveblockingfactor )
|
||||||
|
{
|
||||||
|
throw std::runtime_error("nyx output plug-in requires that \'blocking_factor\' is set!");
|
||||||
|
}
|
||||||
|
|
||||||
|
the_sim_header.dimensions.push_back( 1<<levelmin_ );
|
||||||
|
the_sim_header.dimensions.push_back( 1<<levelmin_ );
|
||||||
|
the_sim_header.dimensions.push_back( 1<<levelmin_ );
|
||||||
|
|
||||||
|
the_sim_header.offset.push_back( 0 );
|
||||||
|
the_sim_header.offset.push_back( 0 );
|
||||||
|
the_sim_header.offset.push_back( 0 );
|
||||||
|
|
||||||
|
the_sim_header.a_start = 1.0/(1.0+cf.get_value<double>("setup","zstart"));
|
||||||
|
the_sim_header.dx = cf.get_value<double>("setup","boxlength")/the_sim_header.dimensions[0]/(cf.get_value<double>("cosmology","H0")*0.01); // not sure?!?
|
||||||
|
the_sim_header.boxlength=cf.get_value<double>("setup","boxlength");
|
||||||
|
the_sim_header.h0 = cf.get_value<double>("cosmology","H0")*0.01;
|
||||||
|
|
||||||
|
if( bhave_hydro )
|
||||||
|
the_sim_header.omega_b = cf.get_value<double>("cosmology","Omega_b");
|
||||||
|
else
|
||||||
|
the_sim_header.omega_b = 0.0;
|
||||||
|
|
||||||
|
the_sim_header.omega_m = cf.get_value<double>("cosmology","Omega_m");
|
||||||
|
the_sim_header.omega_v = cf.get_value<double>("cosmology","Omega_L");
|
||||||
|
the_sim_header.vfact = cf.get_value<double>("cosmology","vfact")*the_sim_header.h0; //.. need to multiply by h, nyx wants this factor for non h-1 units
|
||||||
|
|
||||||
|
std::cout << "creating output object" << std::endl;
|
||||||
|
|
||||||
|
|
||||||
|
//Fix these!
|
||||||
|
lunit_ = the_sim_header.boxlength;
|
||||||
|
vunit_ = the_sim_header.vfact;
|
||||||
|
munit_ = 1.0;
|
||||||
|
|
||||||
|
}
|
||||||
|
|
||||||
|
//destructor
|
||||||
|
virtual ~nyx_output_plugin()
|
||||||
|
{
|
||||||
|
std::string FullPath = fname_;
|
||||||
|
if (!amrex::UtilCreateDirectory(FullPath, 0755))
|
||||||
|
amrex::CreateDirectoryFailed(FullPath);
|
||||||
|
if (!FullPath.empty() && FullPath[FullPath.size()-1] != '/')
|
||||||
|
FullPath += '/';
|
||||||
|
FullPath += "Header";
|
||||||
|
std::ofstream Header(FullPath.c_str());
|
||||||
|
|
||||||
|
for(int lev=0; lev <= f_lev; lev++)
|
||||||
|
{
|
||||||
|
writeLevelPlotFile ( fname_,
|
||||||
|
Header,
|
||||||
|
amrex::VisMF::OneFilePerCPU,
|
||||||
|
lev);
|
||||||
|
//FIXME I would prefer VisMF::NFiles
|
||||||
|
}
|
||||||
|
Header.close();
|
||||||
|
|
||||||
|
writeGridsFile(fname_);
|
||||||
|
std::cout << "destroying output object" << std::endl;
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
void write_particle_data(const particle::container &pc, const cosmo_species &s, double Omega_species ) {};
|
||||||
|
|
||||||
|
// void write_grid_data(const Grid_FFT<real_t> &g, const cosmo_species &s, const fluid_component &c);
|
||||||
|
|
||||||
|
output_type write_species_as(const cosmo_species &s) const
|
||||||
|
{
|
||||||
|
if (s == cosmo_species::baryon)
|
||||||
|
return output_type::field_eulerian;
|
||||||
|
return output_type::field_lagrangian;
|
||||||
|
}
|
||||||
|
|
||||||
|
bool has_64bit_reals() const{ return false; }
|
||||||
|
|
||||||
|
bool has_64bit_ids() const{ return false; }
|
||||||
|
|
||||||
|
real_t position_unit() const { return lunit_; }
|
||||||
|
|
||||||
|
real_t velocity_unit() const { return vunit_; }
|
||||||
|
|
||||||
|
real_t mass_unit() const { return munit_; }
|
||||||
|
|
||||||
|
|
||||||
|
void finalize( void )
|
||||||
|
{
|
||||||
|
//
|
||||||
|
//before finalizing we write out an inputs and a probin file for Nyx.
|
||||||
|
//
|
||||||
|
std::ofstream inputs("inputs");
|
||||||
|
std::ofstream probin("probin");
|
||||||
|
|
||||||
|
//at first the fortran stuff...
|
||||||
|
probin << "&fortin" << std::endl;
|
||||||
|
probin << " comoving_OmM = " << the_sim_header.omega_m << "d0" << std::endl;
|
||||||
|
probin << " comoving_OmB = " << the_sim_header.omega_b << "d0" << std::endl;
|
||||||
|
probin << " comoving_OmL = " << the_sim_header.omega_v << "d0" << std::endl;
|
||||||
|
probin << " comoving_h = " << the_sim_header.h0 << "d0" << std::endl;
|
||||||
|
probin << "/" << std::endl;
|
||||||
|
probin << std::endl;
|
||||||
|
|
||||||
|
//afterwards the cpp stuff...(for which we will need a template, which is read in by the code...)
|
||||||
|
inputs << "nyx.final_a = 1.0 " << std::endl;
|
||||||
|
inputs << "max_step = 100000 " << std::endl;
|
||||||
|
inputs << "nyx.small_dens = 1e-4" << std::endl;
|
||||||
|
inputs << "nyx.small_temp = 10" << std::endl;
|
||||||
|
inputs << "nyx.cfl = 0.9 # cfl number for hyperbolic system" << std::endl;
|
||||||
|
inputs << "nyx.init_shrink = 1.0 # scale back initial timestep" << std::endl;
|
||||||
|
inputs << "nyx.change_max = 1.05 # scale back initial timestep" << std::endl;
|
||||||
|
inputs << "nyx.dt_cutoff = 5.e-20 # level 0 timestep below which we halt" << std::endl;
|
||||||
|
inputs << "nyx.sum_interval = 1 # timesteps between computing mass" << std::endl;
|
||||||
|
inputs << "nyx.v = 1 # verbosity in Castro.cpp" << std::endl;
|
||||||
|
inputs << "gravity.v = 1 # verbosity in Gravity.cpp" << std::endl;
|
||||||
|
inputs << "amr.v = 1 # verbosity in Amr.cpp" << std::endl;
|
||||||
|
inputs << "mg.v = 0 # verbosity in Amr.cpp" << std::endl;
|
||||||
|
inputs << "particles.v = 1 # verbosity in Particle class" << std::endl;
|
||||||
|
inputs << "amr.ref_ratio = 2 2 2 2 2 2 2 2 " << std::endl;
|
||||||
|
inputs << "amr.regrid_int = 2 2 2 2 2 2 2 2 " << std::endl;
|
||||||
|
inputs << "amr.initial_grid_file = init/grids_file" << std::endl;
|
||||||
|
inputs << "amr.useFixedCoarseGrids = 1" << std::endl;
|
||||||
|
inputs << "amr.check_file = chk " << std::endl;
|
||||||
|
inputs << "amr.check_int = 10 " << std::endl;
|
||||||
|
inputs << "amr.plot_file = plt " << std::endl;
|
||||||
|
inputs << "amr.plot_int = 10 " << std::endl;
|
||||||
|
inputs << "amr.derive_plot_vars = particle_count particle_mass_density pressure" << std::endl;
|
||||||
|
inputs << "amr.plot_vars = ALL" << std::endl;
|
||||||
|
inputs << "nyx.add_ext_src = 0" << std::endl;
|
||||||
|
inputs << "gravity.gravity_type = PoissonGrav " << std::endl;
|
||||||
|
inputs << "gravity.no_sync = 1 " << std::endl;
|
||||||
|
inputs << "gravity.no_composite = 1 " << std::endl;
|
||||||
|
inputs << "mg.bottom_solver = 1 " << std::endl;
|
||||||
|
inputs << "geometry.is_periodic = 1 1 1 " << std::endl;
|
||||||
|
inputs << "geometry.coord_sys = 0 " << std::endl;
|
||||||
|
inputs << "amr.max_grid_size = 32 " << std::endl;
|
||||||
|
inputs << "nyx.lo_bc = 0 0 0 " << std::endl;
|
||||||
|
inputs << "nyx.hi_bc = 0 0 0 " << std::endl;
|
||||||
|
inputs << "nyx.do_grav = 1 " << std::endl;
|
||||||
|
inputs << "nyx.do_dm_particles = 1 " << std::endl;
|
||||||
|
inputs << "nyx.particle_init_type = Cosmological " << std::endl;
|
||||||
|
inputs << "nyx.print_fortran_warnings = 0" << std::endl;
|
||||||
|
inputs << "cosmo.initDirName = init " << std::endl;
|
||||||
|
inputs << "nyx.particle_move_type = Gravitational" << std::endl;
|
||||||
|
inputs << "amr.probin_file = probin " << std::endl;
|
||||||
|
inputs << "cosmo.ic-source = MUSIC " << std::endl;
|
||||||
|
|
||||||
|
|
||||||
|
inputs << "amr.blocking_factor = " << cf_.get_value<double>("setup","blocking_factor") << std::endl;
|
||||||
|
|
||||||
|
inputs << "nyx.do_hydro = "<< (the_sim_header.omega_b>0?1:0) << std::endl;
|
||||||
|
inputs << "amr.max_level = " << levelmax_-levelmin_ << std::endl;
|
||||||
|
inputs << "nyx.initial_z = " << 1/the_sim_header.a_start-1 << std::endl;
|
||||||
|
inputs << "amr.n_cell = " << cf_.get_value<int>("setup", "GridRes") << " " << cf_.get_value<int>("setup", "GridRes") << " " << cf_.get_value<int>("setup", "GridRes") << std::endl;
|
||||||
|
inputs << "nyx.n_particles = " << cf_.get_value<int>("setup", "GridRes") << " " << cf_.get_value<int>("setup", "GridRes") << " " << cf_.get_value<int>("setup", "GridRes") << std::endl;
|
||||||
|
inputs << "geometry.prob_lo = 0 0 0" << std::endl;
|
||||||
|
|
||||||
|
//double dx = the_sim_header.dx/the_sim_header.h0;
|
||||||
|
double bl = the_sim_header.boxlength/the_sim_header.h0;
|
||||||
|
inputs << "geometry.prob_hi = " << bl << " " << bl << " " << bl << std::endl;
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
probin.close();
|
||||||
|
inputs.close();
|
||||||
|
std::cout << "finalizing..." << std::endl;
|
||||||
|
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
void writeLevelPlotFile (const std::string& dir,
|
||||||
|
std::ostream& os,
|
||||||
|
amrex::VisMF::How how,
|
||||||
|
int level)
|
||||||
|
{
|
||||||
|
int i, n;
|
||||||
|
|
||||||
|
std::cout << "in writeLevelPlotFile" << std::endl;
|
||||||
|
double h0 = cf_.get_value<double>("cosmology", "H0")*0.01;
|
||||||
|
|
||||||
|
|
||||||
|
if (level == 0)
|
||||||
|
{
|
||||||
|
//
|
||||||
|
// The first thing we write out is the plotfile type.
|
||||||
|
//
|
||||||
|
os << "MUSIC_for_Nyx_v0.1" << '\n';
|
||||||
|
|
||||||
|
os << n_data_items << '\n';
|
||||||
|
|
||||||
|
for (i = 0; i < n_data_items; i++)
|
||||||
|
os << field_name[i] << '\n';
|
||||||
|
|
||||||
|
os << 3 << '\n';
|
||||||
|
os << 0 << '\n';
|
||||||
|
|
||||||
|
os << f_lev << '\n';
|
||||||
|
|
||||||
|
for (i = 0; i < BL_SPACEDIM; i++)
|
||||||
|
os << 0 << ' '; //ProbLo
|
||||||
|
os << '\n';
|
||||||
|
double boxlength = cf_.get_value<double>("setup","boxlength");
|
||||||
|
for (i = 0; i < BL_SPACEDIM; i++)
|
||||||
|
os << boxlength/h0 << ' '; //ProbHi
|
||||||
|
os << '\n';
|
||||||
|
|
||||||
|
for (i = 0; i < f_lev; i++)
|
||||||
|
os << 2 << ' '; //refinement factor
|
||||||
|
os << '\n';
|
||||||
|
|
||||||
|
amrex::IntVect pdLo(0,0,0);
|
||||||
|
amrex::IntVect pdHi(gridp-1,gridp-1,gridp-1);
|
||||||
|
for (i = 0; i <= f_lev; i++) //Geom(i).Domain()
|
||||||
|
{
|
||||||
|
amrex::Box probDomain(pdLo,pdHi);
|
||||||
|
os << probDomain << ' ';
|
||||||
|
pdHi *= 2;
|
||||||
|
pdHi += 1;
|
||||||
|
}
|
||||||
|
os << '\n';
|
||||||
|
|
||||||
|
for (i = 0; i <= f_lev; i++) //level steps
|
||||||
|
os << 0 << ' ';
|
||||||
|
os << '\n';
|
||||||
|
|
||||||
|
double dx = cf_.get_value<double>("setup","boxlength")/gridp/h0;
|
||||||
|
for (i = 0; i <= f_lev; i++)
|
||||||
|
{
|
||||||
|
for (int k = 0; k < BL_SPACEDIM; k++)
|
||||||
|
os << dx << ' ';
|
||||||
|
os << '\n';
|
||||||
|
dx = dx/2.;
|
||||||
|
}
|
||||||
|
os << 0 << '\n';
|
||||||
|
os << "0\n"; // Write bndry data.
|
||||||
|
}
|
||||||
|
|
||||||
|
//
|
||||||
|
// Build the directory to hold the MultiFab at this level.
|
||||||
|
// The name is relative to the directory containing the Header file.
|
||||||
|
//
|
||||||
|
static const std::string BaseName = "/Cell";
|
||||||
|
|
||||||
|
std::string Level = amrex::Concatenate("Level_", level, 1);
|
||||||
|
//
|
||||||
|
// Now for the full pathname of that directory.
|
||||||
|
//
|
||||||
|
std::string FullPath = dir;
|
||||||
|
if (!FullPath.empty() && FullPath[FullPath.size()-1] != '/')
|
||||||
|
FullPath += '/';
|
||||||
|
FullPath += Level;
|
||||||
|
//
|
||||||
|
// Only the I/O processor makes the directory if it doesn't already exist.
|
||||||
|
//
|
||||||
|
if (!amrex::UtilCreateDirectory(FullPath, 0755))
|
||||||
|
amrex::CreateDirectoryFailed(FullPath);
|
||||||
|
|
||||||
|
os << level << ' ' << boxarrays[level].size() << ' ' << 0 << '\n';
|
||||||
|
os << 0 << '\n';
|
||||||
|
|
||||||
|
double cellsize[3];
|
||||||
|
double dx = cf_.get_value<double>("setup","boxlength")/gridp/h0;
|
||||||
|
for (n = 0; n < BL_SPACEDIM; n++)
|
||||||
|
{
|
||||||
|
cellsize[n] = dx;
|
||||||
|
}
|
||||||
|
for (i = 0; i < level; i++)
|
||||||
|
{
|
||||||
|
for (n = 0; n < BL_SPACEDIM; n++)
|
||||||
|
{
|
||||||
|
cellsize[n] /= 2.;
|
||||||
|
}
|
||||||
|
}
|
||||||
|
std::cout << cellsize[0] << std::endl;
|
||||||
|
for (i = 0; i < boxarrays[level].size(); ++i)
|
||||||
|
{
|
||||||
|
double problo[] = {0,0,0};
|
||||||
|
std::cout << boxarrays[level][i] << std::endl;
|
||||||
|
amrex::RealBox gridloc = amrex::RealBox(boxarrays[level][i], cellsize, problo);
|
||||||
|
for (n = 0; n < BL_SPACEDIM; n++)
|
||||||
|
os << gridloc.lo(n) << ' ' << gridloc.hi(n) << '\n';
|
||||||
|
}
|
||||||
|
//
|
||||||
|
// The full relative pathname of the MultiFabs at this level.
|
||||||
|
// The name is relative to the Header file containing this name.
|
||||||
|
// It's the name that gets written into the Header.
|
||||||
|
//
|
||||||
|
std::string PathNameInHeader = Level;
|
||||||
|
PathNameInHeader += BaseName;
|
||||||
|
os << PathNameInHeader << '\n';
|
||||||
|
|
||||||
|
//
|
||||||
|
// Use the Full pathname when naming the MultiFab.
|
||||||
|
//
|
||||||
|
std::string TheFullPath = FullPath;
|
||||||
|
TheFullPath += BaseName;
|
||||||
|
amrex::VisMF::Write((*mfs[level]),TheFullPath,how,true);
|
||||||
|
}
|
||||||
|
|
||||||
|
void writeGridsFile (const std::string& dir)
|
||||||
|
{
|
||||||
|
int i;
|
||||||
|
|
||||||
|
std::cout << "in writeGridsFile" << std::endl;
|
||||||
|
|
||||||
|
std::string myFname = dir;
|
||||||
|
if (!myFname.empty() && myFname[myFname.size()-1] != '/')
|
||||||
|
myFname += '/';
|
||||||
|
myFname += "grids_file";
|
||||||
|
|
||||||
|
std::ofstream os(myFname.c_str());
|
||||||
|
|
||||||
|
os << f_lev << '\n';
|
||||||
|
|
||||||
|
for (int lev = 1; lev <= f_lev; lev++)
|
||||||
|
{
|
||||||
|
os << boxarrays[lev].size() << '\n';
|
||||||
|
boxarrays[lev].coarsen(2);
|
||||||
|
for (i=0; i < boxarrays[lev].size(); i++)
|
||||||
|
os << boxarrays[lev][i] << "\n";
|
||||||
|
}
|
||||||
|
os.close();
|
||||||
|
}
|
||||||
|
};
|
||||||
|
namespace{
|
||||||
|
output_plugin_creator_concrete<nyx_output_plugin> creator("nyx");
|
||||||
|
}
|
||||||
|
|
||||||
|
#endif //ENABLE_AMREX
|
||||||
|
|
Loading…
Reference in a new issue