Skip to content

AMRVAC dat file v5 compatibility; bug-fix on Morton encoding; other m… #26

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
krishnadotprasad wants to merge 1 commit into
base: master
Choose a base branch
from
Open

AMRVAC dat file v5 compatibility; bug-fix on Morton encoding; other m… #26

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion fomo-c/docfiles/amrvac.dox
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -35,7 +35,7 @@ flexibility and removes the need for compilation after each change of option.
For rendering the AMRVAC files with the AIA 193 filter and putting the results in a separate directory /path/to/datfiles/aia193, I would run the following
command
\code{.sh}
example/example_mpi_amrvac/render_all_datfiles -f /path/to/datfiles/file\*.dat -p /path/to/datfiles/amrvac.par -n 1e9 -T 1e6 -L 1e1 -o /path/to/datfiles/aia193 -c ../chiantitables/goft_tables_aia193_abco.dat -x 400 -y 400 -z 600 --lambda_pixel 1 -l 0.524 -b 1.047
example/example_mpi_amrvac/render_all_datfiles -f /path/to/datfiles/file\*.dat -p /path/to/datfiles/amrvac.par -n 1e9 -T 1e6 -L 1e1 -o /path/to/datfiles/aia193 -c ../chiantitables/goft_table_aia193_abco.dat -x 400 -y 400 -z 600 --lambda_pixel 1 -l 0.524 -b 1.047
\endcode
Let me explain these options separately. The first option
\code{.sh}
Expand Down
51 changes: 42 additions & 9 deletions fomo-c/example/example_mpi_amrvac/read_amrvac_files.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -15,7 +15,9 @@ double kboltz = GSL_CONST_MKSA_BOLTZMANN; // # J/K = kg/K m_^2
double mhydro = GSL_CONST_MKSA_MASS_PROTON; // # kg (proton mass for the amrvac_nonAMR.par)
double R_spec = kboltz/(0.5e0*mhydro);

std::vector<int> compatible_versions={3};
//Version 5 (AMRVAC 2.2) compatibility added by Krishna
//Version 4 may also be integrated here after testing
std::vector<int> compatible_versions={3, 5};

// this recursive procedure loops through the forest and calculates the indices of each leaf block
// the implementation is a carbon copy of the python code of Oliver Porth (included in AMRVAC)
Expand Down Expand Up @@ -175,6 +177,16 @@ std::vector<std::vector<double>> read_blocks(istream & file, const int ndim, con
}
}
///////////////////////////////////////////////////////////////
// added by Krishna for compatibility with AMRVAC 2.2 (version 5)
///////////////////////////////////////////////////////////////
if (version==5)
{
for (int j=0; j<ndim; j++)
{
file.read(reinterpret_cast<char*>(&tempint), sizeof(int));
file.read(reinterpret_cast<char*>(&tempint), sizeof(int));
}
}
///////////////////////////////////////////////////////////////
// then start reading the block
// in principle, if the above number of ghost cells is non-zero, then we need to take a different block size, and select different numbers
Expand Down Expand Up @@ -205,7 +217,11 @@ std::vector<std::vector<int>> build_block_info_morton(std::vector<int> nblocks,
// later on, we will prune simulation points from the mortoncurve
// can this be optimized for 2D data?
// If number of blocks is odd, make cube with +1 size
if (maxgridlength % 2 != 0) maxgridlength++;
//if (maxgridlength % 2 != 0) maxgridlength++; //removed by Krishna
//it seems mortonEncode_for will produce indices of the same size as the number of grid cells only
//when the dimensions of the grid is in powers of 2. So added the below line to accommodate this.
//perhaps not an efficient approach when maxgridlength is large. Should be optimized in the future.
maxgridlength = pow(2, ceil(log(maxgridlength)/log(2))); //added by Krishna
mortoncurve.resize(pow(maxgridlength,3));
for (int k=0; k<maxgridlength; k++)
for (int j=0; j<maxgridlength; j++)
Expand Down Expand Up @@ -488,8 +504,26 @@ FoMo::FoMoObject read_amrvac_new_dat_file(const char* datfile, int version_numbe
file.read(reinterpret_cast<char*>(&nx.at(i)), sizeof(int));
nglev1*=nx.at(i);
std::cout << "block_nx(" << i << "):" << nx.at(i) << std::endl;
}
///////////////////////////////////////////////////////////////
// added by Krishna for compatibility with AMRVAC 2.2 (version 5)
///////////////////////////////////////////////////////////////
if (version == 5)
{
std::vector<int> periodic(ndim);
for (int i=0; i<ndim; i++)
{
file.read(reinterpret_cast<char*>(&periodic.at(i)), sizeof(int));
std::cout << "Boundary Periodicity ["<< i <<"]:" << periodic.at(i) << std::endl;
}
char geom[16];
file.read(geom,16);
std::cout << "Geometry:" << geom <<std::endl;
int stag;
file.read(reinterpret_cast<char*>(&stag),sizeof(int));
std::cout << "Staggered:" << stag <<std::endl;
}

///////////////////////////////////////////////////////////////
char temp[16];
std::vector<std::string> varnames(nw);
for (int i=0; i<nw; i++)
Expand All @@ -510,7 +544,7 @@ FoMo::FoMoObject read_amrvac_new_dat_file(const char* datfile, int version_numbe
double gamma=0;
std::vector<double> eqpar(neqpar);
std::vector<std::string> eqparnames(neqpar);
std::string gammaname="gamma ";
std::string gammaname="gamma"; //removed spaces by Krishna for compatibility with versions 3 & 5
for (int i=0; i<neqpar; i++)
{
file.read(reinterpret_cast<char*>(&tmpeqpar), sizeof(double));
Expand All @@ -522,7 +556,8 @@ FoMo::FoMoObject read_amrvac_new_dat_file(const char* datfile, int version_numbe
eqparnames.at(i)=temp;
// if the variable name matches gammaname, then take that value for gamma
// this has as consequence that the last of such values will eventually be used
if (eqparnames.at(i).compare(gammaname)==0) gamma=eqpar.at(i);
if (eqparnames.at(i).substr(0,5).compare(gammaname)==0) gamma=eqpar.at(i); //added substr by Krishna
std::cout << "Gamma:" << gamma << std::endl; //added by Krishna
}
if (gamma==0.)
{
Expand All @@ -544,6 +579,7 @@ FoMo::FoMoObject read_amrvac_new_dat_file(const char* datfile, int version_numbe
//compute the dx in each direction
std::vector<double> cellsize(ndim);
blocksandsize(nxlone,nx,xprobmin,xprobmax,nblocks,cellsize);
for (int i=0; i<ndim; i++) std::cout << "cell size ["<< i << "]:" << cellsize[i] << std::endl; //added by Krishna

// make sure the code is also working for 2D data!
if (ndim<3) nblocks.push_back(1);
Expand Down Expand Up @@ -728,8 +764,6 @@ FoMo::FoMoObject read_amrvac_dat_file_v4(const char* datfile, int version_number
//////////////////////////////////////////////////////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////////////////////////////////////////////////////



FoMo::FoMoObject read_amrvac_dat_file(const char* datfile, const char* amrvacpar, string amrvac_version, const int gamma_eqparposition, const double n_unit, const double Teunit, const double L_unit)
{
int version_number=read_dat_version(datfile);
Expand All @@ -742,13 +776,12 @@ FoMo::FoMoObject read_amrvac_dat_file(const char* datfile, const char* amrvacpar
Object=read_amrvac_new_dat_file(datfile,version_number,n_unit,Teunit,L_unit);
}
//////////// Added by vaibhav pant for DATfiles version 4///////////////////////////////////////////////////
if(version_number == 4) {
else if(version_number == 4) { //modified (to "else if") by Krishna
/* We use the new snapshot reader, which is similar to version 3 except few changes in the definition of Gamma */
std::cout << "Using new datfile reader (version 4) for AMRVAC 2.0" << std::endl;
Object=read_amrvac_dat_file_v4(datfile,version_number,n_unit,Teunit,L_unit);
}
////////////////////////////////////////////////////////////////////////////////////////////////////////////

else {
/* Not a compatible version, use the old reader */
std::cout << "Using old datfile reader." << std::endl;
Expand Down
16 changes: 8 additions & 8 deletions fomo-c/python/amrvac_animation.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,24 +2,23 @@
import numpy as np
import pylab as plt
import matplotlib.cm as cm
from readfomo import readgoftcube, readgoftcubechianti, regulargoftcube, gaussfitgoftcube
from readfomo import readgoftcube, readgoftcube_units_chianti, regulargoftcube, gaussfitgoftcube
import glob
import sys
import os
import sunpy.cm

import sunpy.visualization.colormaps.cm

# store filenames as a sorted list
unsorted=glob.glob(sys.argv[1])
filelist=sorted(unsorted)

# create data cube from file list
print 'Reading file',filelist[0]
data,chiantifile=readgoftcubechianti(filelist[0])
print("Reading file:",filelist[0])
data,units,chiantifile=readgoftcube_units_chianti(filelist[0]) #modified by Krishna
emiss,xvec,yvec,lvec=regulargoftcube(data)
allemiss=emiss
for filename in filelist[1:]:
print 'Reading file',filename
print("Reading file:",filename)
data=readgoftcube(filename)
emiss,xvec,yvec,lvec=regulargoftcube(data)
allemiss=np.dstack((allemiss,emiss))
Expand All @@ -28,12 +27,13 @@
fig = plt.figure(figsize=(10,10))
ax = fig.add_subplot(111)
ax.set_aspect('auto')
ax.set_xlabel('x (Mm)')
ax.set_ylabel('y (Mm)')
ax.set_xlabel('x ('+units[0]+')')
ax.set_ylabel('y ('+units[1]+')')

# set colour map
cmap=cm.hot
aiafilter=os.path.basename(chiantifile)
aiafilter=aiafilter.decode('ascii') #added by Krishna
if ('goft_table_aia094' in aiafilter): cmap=plt.get_cmap('sdoaia094')
if ('goft_table_aia131' in aiafilter): cmap=plt.get_cmap('sdoaia131')
if ('goft_table_aia171' in aiafilter): cmap=plt.get_cmap('sdoaia171')
Expand Down