Skip to content
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

Remove VTK_cartesian config key and corresponding code #34

Merged
merged 2 commits into from
Jan 6, 2025
Merged

Remove VTK_cartesian config key and corresponding code #34

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
df6e7d3
Remove code duplications in VTK source code
RobinSattler Dec 30, 2024
b8c0fb0
Remove VTK_cartesian key and corresponding code
RobinSattler Jan 3, 2025
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
8 changes: 3 additions & 5 deletions src/main.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -69,7 +69,7 @@ void checkGridBorders(double min, double max, std::string _x) {

int nx {100}, ny {100}, nz {100}, eosType {1}, etaSparam {0}, zetaSparam {0}, eosTypeHadron {0};
// only FO hypersurface output: {0,1}; freezeout output extended by e,nb: {0,1}
bool vtk_cartesian {false}, freezeoutOnly {false}, freezeoutExtend {false}, vorticityOn {false};
bool freezeoutOnly {false}, freezeoutExtend {false}, vorticityOn {false};
double xmin {-5.0}, xmax {5.0}, ymin {-5.0}, ymax {5.0}, etamin {-5.0},
etamax {5.0}, tau0 {1.0}, tauMax {20.0}, tauResize {4.0}, dtau {0.05},
etaS {0.08}, zetaS {0.0}, eCrit {0.5}, etaSEpsilonMin {5.}, al {0.}, ah {0.}, aRho {0.}, T0 {0.15},
Expand Down Expand Up @@ -121,7 +121,6 @@ void readParameters(char *parFile) {
{"impactPar", [](const string& value) { impactPar = atof(value.c_str()); }},
{"s0ScaleFactor", [](const string& value) { s0ScaleFactor = atof(value.c_str()); }},
{"VTK_output_values", [](const string& value) { vtk_values = value; }},
{"VTK_cartesian", [](const string& value) { vtk_cartesian= atoi(value.c_str()); }},
{"etaSparam", [](const string& value) { etaSparam = atoi(value.c_str()); }},
{"aRho", [](const string& value) { aRho = atof(value.c_str()); }},
{"ah", [](const string& value) { ah = atof(value.c_str()); }},
Expand Down Expand Up @@ -220,8 +219,7 @@ void printParameters() {
cout << "smoothingType = " << smoothingType << endl;
cout << "impactPar = " << impactPar << endl;
cout << "s0ScaleFactor = " << s0ScaleFactor << endl;
cout << "VTK output values = " << vtk_values << endl;
cout << "VTK cartesian = " << vtk_cartesian << endl;
cout << "VTK_output_values = " << vtk_values << endl;
cout << "======= end parameters =======\n";
}

Expand Down Expand Up @@ -395,7 +393,7 @@ int main(int argc, char **argv) {
bool resized = false; // flag if the grid has been resized

std::string dir=outputDir.c_str();
VtkOutput vtk_out=VtkOutput(dir,eos,xmin,ymin,etamin,vtk_cartesian);
VtkOutput vtk_out=VtkOutput(dir,eos,xmin,ymin,etamin);

int nelements = 0;
do {
Expand Down
201 changes: 58 additions & 143 deletions src/vtk.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -8,26 +8,22 @@
#include "hdo.h"
#include "vtk.h"

void VtkOutput::write_header(std::ofstream &file, const Hydro h, const std::string &description) {
if (cartesian_) {
file << "# vtk DataFile Version 2.0\n"
<< description << "\n"
<< "ASCII\n"
<< "DATASET STRUCTURED_POINTS\n"
<< "DIMENSIONS " << h.getFluid()->getNX() << " " << h.getFluid()->getNY() << " " << (h.getFluid()->getNZ())*20 << "\n"
<< "SPACING " << h.getFluid()->getDx() << " " << h.getFluid()->getDy() << " " << h.getTau()*std::sinh(h.getFluid()->getDz()*(h.getFluid()->getNZ()/2))/(20*h.getFluid()->getNZ()/2) << "\n"
<< "ORIGIN " << xmin_ << " " << ymin_ << " " << etamin_ << "\n"
<< "POINT_DATA " << h.getFluid()->getNX()*h.getFluid()->getNY()*(h.getFluid()->getNZ())*20 << "\n";
} else {
file << "# vtk DataFile Version 2.0\n"
<< description << "\n"
<< "ASCII\n"
<< "DATASET STRUCTURED_POINTS\n"
<< "DIMENSIONS " << h.getFluid()->getNX() << " " << h.getFluid()->getNY() << " " << h.getFluid()->getNZ() << "\n"
<< "SPACING " << h.getFluid()->getDx() << " " << h.getFluid()->getDy() << " " << h.getFluid()->getDz() << "\n"
<< "ORIGIN " << xmin_ << " " << ymin_ << " " << etamin_ << "\n"
<< "POINT_DATA " << h.getFluid()->getNX()*h.getFluid()->getNY()*h.getFluid()->getNZ() << "\n";
}
void VtkOutput::write_header(std::ofstream &file, const Hydro h,
const std::string &description) {
num_of_cells_x_direction_ = h.getFluid()->getNX();
num_of_cells_y_direction_ = h.getFluid()->getNY();
num_of_cells_eta_direction_ = h.getFluid()->getNZ();
file << "# vtk DataFile Version 2.0\n"
<< description << "\n"
<< "ASCII\n"
<< "DATASET STRUCTURED_POINTS\n"
<< "DIMENSIONS " << num_of_cells_x_direction_ << " "
<< num_of_cells_y_direction_ << " " << num_of_cells_eta_direction_ << "\n"
<< "SPACING " << h.getFluid()->getDx() << " " << h.getFluid()->getDy()
<< " " << h.getFluid()->getDz() << "\n"
<< "ORIGIN " << xmin_ << " " << ymin_ << " " << etamin_ << "\n"
<< "POINT_DATA " << num_of_cells_x_direction_ * num_of_cells_y_direction_
* num_of_cells_eta_direction_ << "\n";

return;
}
Expand All @@ -38,98 +34,47 @@ std::string VtkOutput::make_filename(const std::string &descr, int counter) {
return path_ + std::string("/") + descr + std::string(suffix);
}

std::vector<double> VtkOutput::smearing_factor_and_poseta(const Hydro h, const int iz, const int z_length) {
double poseta = iz;
double factor = 1;
if (cartesian_) {
double total_length = h.getTau()*std::sinh(h.getFluid()->getDz()*(h.getFluid()->getNZ()/2));
double pos = 0;
if (iz < z_length/2) {
pos = -total_length+total_length/(z_length/2)*iz;
for (int ieta = 0; ieta < h.getFluid()->getNZ()/2; ieta++) {
if (h.getTau()*std::sinh(h.getFluid()->getDz()*(ieta-h.getFluid()->getNZ()/2)) > pos) {
factor = fabs((h.getTau()*std::sinh(h.getFluid()->getDz()*(ieta-h.getFluid()->getNZ()/2))-pos)/pos);
poseta = ieta;
break;
}
}
} else {
pos = total_length/(z_length/2)*(iz-z_length/2);
for (int ieta = h.getFluid()->getNZ()/2; ieta < h.getFluid()->getNZ(); ieta++) {
if (h.getTau()*std::sinh(h.getFluid()->getDz()*(ieta-h.getFluid()->getNZ()/2)) > pos) {
factor = fabs((h.getTau()*std::sinh(h.getFluid()->getDz()*(ieta-h.getFluid()->getNZ()/2))-pos)/pos);
poseta = ieta;
break;
}
}
}
if (pos == 0) {
factor = 1;
}
}
std::vector<double> factor_and_poseta = {factor, poseta};
return factor_and_poseta;
}

void VtkOutput::write_vtk_scalar(std::ofstream &file, const Hydro h,
std::string &quantity) {
const std::string &quantity) {
file << "SCALARS " << quantity << " double 1\n"
<< "LOOKUP_TABLE default\n";
file << std::setprecision(3);
file << std::fixed;
int z_length = h.getFluid()->getNZ();
if (cartesian_) {
z_length = 20*z_length;
}

for (int iz = 0; iz < z_length; iz++) {
std::vector<double> factor_and_poseta = smearing_factor_and_poseta(h, iz, z_length);
double factor = factor_and_poseta.at(0);
double poseta = factor_and_poseta.at(1);
for (int iy = 0; iy < h.getFluid()->getNY(); iy++) {
for (int ix = 0; ix < h.getFluid()->getNX(); ix++) {
double e, mub, muq, mus, nb, nq, ns, p, T, vx, vy, vz;
double e2, mub2, muq2, mus2, nb2, nq2, ns2, p2, T2, vx2, vy2, vz2;
Cell* cell = h.getFluid()->getCell(ix,iy,poseta);
Cell* cell2;
if (poseta > 0) {
cell2 = h.getFluid()->getCell(ix,iy,poseta-1);
} else {
cell2 = cell;
}
for (int ieta = 0; ieta < num_of_cells_eta_direction_; ieta++) {
for (int iy = 0; iy < num_of_cells_y_direction_; iy++) {
for (int ix = 0; ix < num_of_cells_x_direction_; ix++) {
double e, nb, nq, ns, p, vx, vy, vz;
Cell* cell = h.getFluid()->getCell(ix,iy,ieta);
cell->getPrimVar(eos_, h.getTau(), e, p, nb, nq, ns, vx, vy, vz);
cell2->getPrimVar(eos_, h.getTau(), e2, p2, nb2, nq2, ns2, vx2, vy2, vz2);
double q = 0;
// scalar quantities
if (quantity == "eps") {
q = factor*e+(factor-1)*e2;
} else if (quantity == "mub") {
eos_->eos(e, nb, nq, ns, T, mub, muq, mus, p);
eos_->eos(e2, nb2, nq2, ns2, T2, mub2, muq2, mus2, p2);
q = factor*mub+(factor-1)*mub2;
} else if (quantity == "muq") {
eos_->eos(e, nb, nq, ns, T, mub, muq, mus, p);
eos_->eos(e2, nb2, nq2, ns2, T2, mub2, muq2, mus2, p2);
q = factor*muq+(factor-1)*muq2;
} else if (quantity == "mus") {
eos_->eos(e, nb, nq, ns, T, mub, muq, mus, p);
eos_->eos(e2, nb2, nq2, ns2, T2, mub2, muq2, mus2, p2);
q = factor*mus+(factor-1)*mus2;
q = e;
} else if (quantity == "nb") {
q = factor*nb+(factor-1)*nb2;
q = nb;
} else if (quantity == "nq") {
q = factor*nq+(factor-1)*nq2;
q = nq;
} else if (quantity == "ns") {
q = factor*ns+(factor-1)*ns2;
q = ns;
} else if (quantity == "p") {
q = factor*p+(factor-1)*p2;
q = p;
} else if (quantity == "Pi") {
double Pi = cell->getPi();
double Pi2 = cell2->getPi();
q = factor*Pi+(factor-1)*Pi2;
} else if (quantity == "T") {
q = cell->getPi();
// scalar quantities that need eos()
} else if (quantity == "mub" || quantity == "muq" || quantity == "mus"
|| quantity == "T") {
double mub, muq, mus, T;
eos_->eos(e, nb, nq, ns, T, mub, muq, mus, p);
eos_->eos(e2, nb2, nq2, ns2, T2, mub2, muq2, mus2, p2);
q = factor*T+(factor-1)*T2;
if (quantity == "mub") {
q = mub;
} else if (quantity == "muq") {
q = muq;
} else if (quantity == "mus") {
q = mus;
} else if (quantity == "T") {
q = T;
}
}
file << q << " ";
}
Expand All @@ -139,35 +84,20 @@ void VtkOutput::write_vtk_scalar(std::ofstream &file, const Hydro h,
}

void VtkOutput::write_vtk_vector(std::ofstream &file, const Hydro h,
std::string &quantity) {
const std::string &quantity) {
file << "VECTORS " << quantity << " double\n";
file << std::setprecision(3);
file << std::fixed;
int z_length = h.getFluid()->getNZ();
if (cartesian_) {
z_length = 20*z_length;
}

for (int iz = 0; iz < z_length; iz++) {
std::vector<double> factor_and_poseta = smearing_factor_and_poseta(h, iz, z_length);
double factor = factor_and_poseta.at(0);
double poseta = factor_and_poseta.at(1);
for (int iy = 0; iy < h.getFluid()->getNY(); iy++) {
for (int ix = 0; ix < h.getFluid()->getNX(); ix++) {
double e, p, nb, nq, ns, vx, vy, vz, T, mub, muq, mus;
double e2, p2, nb2, nq2, ns2, vx2, vy2, vz2, T2, mub2, muq2, mus2;
Cell* cell = h.getFluid()->getCell(ix,iy,poseta);
Cell* cell2;
if (poseta > 0) {
cell2 = h.getFluid()->getCell(ix,iy,poseta-1);
} else {
cell2 = cell;
}
for (int ieta = 0; ieta < num_of_cells_eta_direction_; ieta++) {
for (int iy = 0; iy < num_of_cells_y_direction_; iy++) {
for (int ix = 0; ix < num_of_cells_x_direction_; ix++) {
double e, p, nb, nq, ns, vx, vy, vz;
Cell* cell = h.getFluid()->getCell(ix,iy,ieta);
cell->getPrimVar(eos_, h.getTau(), e, p, nb, nq, ns, vx, vy, vz);
cell2->getPrimVar(eos_, h.getTau(), e2, p2, nb2, nq2, ns2, vx2, vy2, vz2);
std::vector<double> q = {0.,0.,0.};
if (quantity == "v") {
q = {factor*vx+(factor-1)*vx2, factor*vy+(factor-1)*vy2, factor*vz+(factor-1)*vz2};
q = {vx, vy, vz};
}
file << q.at(0) << " " << q.at(1) << " " << q.at(2) << "\n";
}
Expand All @@ -176,37 +106,22 @@ void VtkOutput::write_vtk_vector(std::ofstream &file, const Hydro h,
}

void VtkOutput::write_vtk_tensor(std::ofstream &file, const Hydro h,
std::string &quantity) {
const std::string &quantity) {
for (int i = 0; i < 4; i++) {
for (int j = 0; j < 4; j++ ) {
file << "SCALARS " << quantity << std::to_string(i) << std::to_string(j)
<< " double 1\n"
<< "LOOKUP_TABLE default\n";
file << std::setprecision(3);
file << std::fixed;
int z_length = h.getFluid()->getNZ();
if (cartesian_) {
z_length = 20*z_length;
}

for (int iz = 0; iz < z_length; iz++) {
std::vector<double> factor_and_poseta = smearing_factor_and_poseta(h, iz, z_length);
double factor = factor_and_poseta.at(0);
double poseta = factor_and_poseta.at(1);
for (int iy = 0; iy < h.getFluid()->getNY(); iy++) {
for (int ix = 0; ix < h.getFluid()->getNX(); ix++) {
Cell* cell = h.getFluid()->getCell(ix,iy,poseta);
Cell* cell2;
if (poseta > 0) {
cell2 = h.getFluid()->getCell(ix,iy,poseta-1);
} else {
cell2 = cell;
}
for (int ieta = 0; ieta < num_of_cells_eta_direction_; ieta++) {
for (int iy = 0; iy < num_of_cells_y_direction_; iy++) {
for (int ix = 0; ix < num_of_cells_x_direction_; ix++) {
Cell* cell = h.getFluid()->getCell(ix,iy,ieta);
double q = 0;
if (quantity == "pi") {
double pi_ij = cell->getpi(i,j);
double pi2_ij = cell2->getpi(i,j);
q = factor*pi_ij+(factor-1)*pi2_ij;
q = cell->getpi(i,j);
}
file << q << " ";
}
Expand All @@ -219,7 +134,7 @@ void VtkOutput::write_vtk_tensor(std::ofstream &file, const Hydro h,

std::vector<std::string> split (const std::string &s, char delim) {
std::vector<std::string> result;
std::stringstream ss (s);
std::stringstream ss(s);
std::string item;

while (getline(ss, item, delim)) {
Expand All @@ -229,13 +144,13 @@ std::vector<std::string> split (const std::string &s, char delim) {
return result;
}

bool VtkOutput::is_quantity_implemented(std::string &quantity) {
bool VtkOutput::is_quantity_implemented(const std::string &quantity) {
bool quantity_is_valid = (valid_quantities_.find(quantity)
!= valid_quantities_.end());
return quantity_is_valid;
}

void VtkOutput::write(const Hydro h, std::string &quantities) {
void VtkOutput::write(const Hydro h, const std::string &quantities) {
std::vector<std::string> quantities_list = split(quantities,',');
for (std::string q : quantities_list){
if (!is_quantity_implemented(q)) {
Expand Down
23 changes: 10 additions & 13 deletions src/vtk.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -18,10 +18,11 @@ class VtkOutput {
std::string path_;
EoS* eos_;
double xmin_, ymin_, etamin_;
bool cartesian_;
int num_of_cells_x_direction_, num_of_cells_y_direction_,
num_of_cells_eta_direction_;
int vtk_output_counter_ = 0; // Number of vtk output in current event
/*
* The following vector contains all currently supported VTK quantities.
* The following map contains all currently supported VTK quantities.
* If multiple quantities are desired, the delimiter in the config file has
* to be a comma without any whitespaces in between, for example:
* `VTK_output_valus eps,mub,nq,T,v,pi`
Expand Down Expand Up @@ -52,29 +53,25 @@ class VtkOutput {
* \param xmin x coordinate of the first cell (center of fluid cell).
* \param ymin y coordinate of the first cell (center of fluid cell).
* \param etamin eta coordinate of the first cell (center of fluid cell).
* \param cartesian Bool that states if Cartesian output is enabled.
*/
VtkOutput(std::string path, EoS* eos, double xmin, double ymin,
double etamin, bool cartesian=false):
double etamin):
path_(path),
eos_(eos),
xmin_(xmin),
ymin_(ymin),
etamin_(etamin),
cartesian_(cartesian)
etamin_(etamin)
{}

void write(const Hydro h, std::string &quantity);
void write(const Hydro h, const std::string &quantities);
void write_header(std::ofstream &file, const Hydro h,
const std::string &description);
void write_vtk_scalar(std::ofstream &file, const Hydro h,
std::string &quantity);
const std::string &quantity);
void write_vtk_vector(std::ofstream &file, const Hydro h,
std::string &quantity);
const std::string &quantity);
void write_vtk_tensor(std::ofstream &file, const Hydro h,
std::string &quantity);
bool is_quantity_implemented(std::string &quantity);
std::vector<double> smearing_factor_and_poseta(const Hydro h, const int iz,
const int z_length);
const std::string &quantity);
bool is_quantity_implemented(const std::string &quantity);
std::string make_filename (const std::string &descr, int counter);
};