Skip to Content.
Sympa Menu

opal - RE: [Opal] OPAL - central region

opal AT lists.psi.ch

Subject: The OPAL Discussion Forum

List archive

RE: [Opal] OPAL - central region


Chronological Thread 
  • From: "Snuverink Jochem (PSI)" <jochem.snuverink AT psi.ch>
  • To: Sébastien Toussaint <sebastien.toussaint AT uclouvain.be>
  • Cc: opal <opal AT lists.psi.ch>
  • Subject: RE: [Opal] OPAL - central region
  • Date: Fri, 24 Jul 2020 15:17:06 +0000
  • Accept-language: en-US, de-CH

Dear Sébastien Toussaint,

 

Thanks for your mail, I would recommend to subscribe to the OPAL mailing list at https://psilists.ethz.ch/sympa/subscribe/opal

 

Attached is a C-program that converts an ASCII output to the h5 format, it can also be found in the repository under tools/BandRF. You might need to adapt it to your use case as it was written for ANSYS output files.

 

Best regards,

Jochem

 

From: Christof Metzger-Kraus <christof.j.kraus AT gmail.com>
Sent: Freitag, 24. Juli 2020 17:01
To: Sébastien Toussaint <sebastien.toussaint AT uclouvain.be>
Cc: Snuverink Jochem (PSI) <jochem.snuverink AT psi.ch>; Frey Matthias (PSI) <matthias.frey AT psi.ch>
Subject: Re: OPAL - central region

 

Dear Mr. Toussaint,

 

I'm not the expert on OPAL-cycl but the two other developers I've added in CC are. I'm sure one of them can answer your question. 

 

Best regards 

Christof Metzger-Kraus 

 

Sébastien Toussaint <sebastien.toussaint AT uclouvain.be> schrieb am Fr., 24. Juli 2020, 16:35:

Dear Dr. Kraus,

 

Sorry if there exists a more suitable way to contact you or the OPAL development team.

 

Last year, after the talk of Prof. Adelmann at CYC2019 (cape town), We started to use OPAL-cycl to simulate orbits in a cyclotron. Our goal would now be to include the RF field to simulate accurately what is happening in the central region. The electrostatic simulations were realized with COMSOL.

 

The BANDRF option of opal requires a .h5 format that I am unable to export from COMSOL. After a lot of unfruitful attempts, I contact you to know if you have a clue to guide us going from a txt file exported from comsol (Basically x,y,z, Ex,Ey,Ez) to a .h5 format ?

 

Any advice as to how we could proceed is welcome.

 

Best regards,

 

Sébastien Toussaint

/*

Purpose: Convert ANSYS E & B-Field data into H5hut (H5block)
format for usage in OPAL.

Usage: ascii2h5block efield.txt hfield.txt ehfout

To visualize use Visit: https://wci.llnl.gov/codes/visit/

Ch. Wang & A. Adelmann, 2011

ToDo: make it more generic

Modification by Chris van Herwaarden and Hui Zhang

The first three rows of a field map that you wish to combine should look
like this:

int1 int2 int3
int1 int2 int3
int1 int2 int3

the integers are the amount of steps (or different values) in x y and z

*/


#include <cmath>
#include <cstdlib>
#include <fstream>
#include <iostream>
#include <string>

#include "H5hut.h"

int main(int argc,char *argv[]) {

if (argc != 4) {
std::cout << "Wrong number of arguments: ascii2h5block efield.txt (or
\"\") hfield.txt (or \"\") ehfout" << std::endl;
//--commlib mpi" << std::endl;
std::exit(1);
}

// // initialize MPI & H5hut
// MPI_Init (&argc, &argv);
// MPI_Comm comm = MPI_COMM_WORLD;
// int comm_size = 1;
// MPI_Comm_size (comm, &comm_size);
// int comm_rank = 0;
// MPI_Comm_rank (comm, &comm_rank);
// H5AbortOnError ();
// H5SetVerbosityLevel (h5_verbosity);

std::string efin(argv[1]);
std::string hfin(argv[2]);
std::string ehfout(argv[3]);
ehfout += std::string(".h5part");

h5_float64_t freq = 72.615 * 1.0e6;

std::cout << "--------------------------------------------------------"
<< std::endl;
std::cout << "Combine " << efin << " and " << hfin << " to " << ehfout <<
std::endl;

std::cout << "Frequency " << freq << " [Hz]" << std::endl;

std::ifstream finE, finH;
finH.open(hfin);
finE.open(efin);
bool efield = false, hfield = false;
if (finE.is_open())
efield = true;
else if (efin.empty() == false) {
std::cout << "E-field \"" << efin << "\" could not be opened" <<
std::endl;
std::exit(1);
}
if (finH.is_open())
hfield = true;
else if (hfin.empty() == false) {
std::cout << "H-field \"" << hfin << "\" could not be opened" <<
std::endl;
std::exit(1);
}

if (efield == false && hfield == false) {
std::cerr << "Neither E-field \"" << efin
<< "\" nor H-field \"" << hfin
<< "\" could be opened" << std::endl;
std::exit(1);
}

int gridPx = 0, gridPy = 0, gridPz = 0;
char temp[256];
/* Header and grid info */
/* Deletes the first two rows in finE and finH to get rid of header info*/
if (efield) {
finE.getline(temp,256);
finE >> gridPx >> gridPy >> gridPz;
finE.getline(temp,256);
finE.getline(temp,256);
}

int HgridPx = 0, HgridPy = 0, HgridPz = 0;
if (hfield) {
finH.getline(temp,256);
finH >> HgridPx >> HgridPy >> HgridPz;
finH.getline(temp,256);
finH.getline(temp,256);
}

int n = gridPx * gridPy * gridPz; /* number of rows in a column */
int m = HgridPx * HgridPy * HgridPz;

int H5gridPx = std::max(gridPx, HgridPx);
int H5gridPy = std::max(gridPy, HgridPy);
int H5gridPz = std::max(gridPz, HgridPz);

std::cout << "H5block grid" << std::endl;
std::cout << "H5gridPx " << H5gridPx << std::endl;
std::cout << "H5gridPy " << H5gridPy << std::endl;
std::cout << "H5gridPz " << H5gridPz << std::endl;

h5_file_t file = H5OpenFile(ehfout.c_str(), H5_O_WRONLY, H5_PROP_DEFAULT);
if (!file) {
std::cerr << "Could not open output file " << ehfout << std::endl;
std::exit(1);
}

H5SetStep(file, 0);
H5Block3dSetView(file,
0, H5gridPx - 1,
0, H5gridPy - 1,
0, H5gridPz - 1);

std::cout << "number Edata " << n << std::endl;
if (efield) {
h5_float64_t* sEx = new h5_float64_t[n]; /* redefines the sE and sH
variables as arrays */
h5_float64_t* sEy = new h5_float64_t[n];
h5_float64_t* sEz = new h5_float64_t[n];

h5_float64_t* FieldstrengthEz = new h5_float64_t[n]; /* redefines the
fieldstrength variables as arrays */
h5_float64_t* FieldstrengthEx = new h5_float64_t[n];
h5_float64_t* FieldstrengthEy = new h5_float64_t[n];

double* Ex = new double[n]; /* redefines the E and H variables as
arrays */
double* Ey = new double[n];
double* Ez = new double[n];

for (int i = 0; i < n; i++) {
finE >> sEx[i] >> sEy[i] >> sEz[i] >> Ex[i] >> Ey[i] >> Ez[i];
}
finE.close();

h5_float64_t stepEx = (sEx[n-1] - sEx[0]) / (gridPx - 1); /*
calculates the stepsizes of the x,y,z of the efield and hfield*/
h5_float64_t stepEy = (sEy[n-1] - sEy[0]) / (gridPy - 1);
h5_float64_t stepEz = (sEz[n-1] - sEz[0]) / (gridPz - 1);

std::cout << "gridPx " << gridPx << " stepEx " << stepEx << std::endl;
std::cout << "gridPy " << gridPy << " stepEy " << stepEy << std::endl;
std::cout << "gridPz " << gridPz << " stepEz " << stepEz << std::endl;

for (int i = 0; i < gridPz; i++) {
for (int j = 0; j < gridPy; j++) {
for (int k = 0; k < gridPx; k++) {
FieldstrengthEx[k + j * gridPx + i * gridPx * gridPy] =
static_cast<h5_float64_t>(Ex[i + j * gridPz + k * gridPz * gridPy]);
FieldstrengthEy[k + j * gridPx + i * gridPx * gridPy] =
static_cast<h5_float64_t>(Ey[i + j * gridPz + k * gridPz * gridPy]);
FieldstrengthEz[k + j * gridPx + i * gridPx * gridPy] =
static_cast<h5_float64_t>(Ez[i + j * gridPz + k * gridPz * gridPy]);
}
}
}

H5Block3dWriteVector3dFieldFloat64 (
file, /*!< IN: file
handle */
"Efield", /*!< IN: name of
dataset to write */
FieldstrengthEx, /*!< IN: X axis
data */
FieldstrengthEy, /*!< IN: Y axis
data */
FieldstrengthEz /*!< IN: Z axis
data */
);
H5Block3dSetFieldSpacing(file, "Efield", stepEx, stepEy, stepEz);
H5Block3dSetFieldOrigin (file, "Efield", sEx[0], sEy[0], sEz[0]);
}

std::cout << "number Bdata " << m << std::endl;

if (hfield) {

h5_float64_t* sHx = new h5_float64_t[m];
h5_float64_t* sHy = new h5_float64_t[m];
h5_float64_t* sHz = new h5_float64_t[m];

h5_float64_t* FieldstrengthHz = new h5_float64_t[m];
h5_float64_t* FieldstrengthHx = new h5_float64_t[m];
h5_float64_t* FieldstrengthHy = new h5_float64_t[m];

double* Hx = new double[m];
double* Hy = new double[m];
double* Hz = new double[m];

for (int i = 0; i < m; i++) {
finH >> sHx[i] >> sHy[i] >> sHz[i] >> Hx[i] >> Hy[i] >> Hz[i];
}
finH.close();

h5_float64_t stepHx = (sHx[m-1] - sHx[0]) / (HgridPx - 1);
h5_float64_t stepHy = (sHy[m-1] - sHy[0]) / (HgridPy - 1);
h5_float64_t stepHz = (sHz[m-1] - sHz[0]) / (HgridPz - 1);

std::cout << "HgridPx " << HgridPx << " stepHx " << stepHx <<
std::endl;
std::cout << "HgridPy " << HgridPy << " stepHy " << stepHy <<
std::endl;
std::cout << "HgridPz " << HgridPz << " stepHz " << stepHz <<
std::endl;

for (int i = 0; i < HgridPz; i++) {
for (int j = 0; j < HgridPy; j++) {
for (int k = 0; k < HgridPx; k++) {
FieldstrengthHx[k + j * HgridPx + i * HgridPx * HgridPy]
= static_cast<h5_float64_t>((Hx[i + j * HgridPz + k * HgridPz * HgridPy]));
FieldstrengthHy[k + j * HgridPx + i * HgridPx * HgridPy]
= static_cast<h5_float64_t>((Hy[i + j * HgridPz + k * HgridPz * HgridPy]));
FieldstrengthHz[k + j * HgridPx + i * HgridPx * HgridPy]
= static_cast<h5_float64_t>((Hz[i + j * HgridPz + k * HgridPz * HgridPy]));
}
}
}
H5Block3dWriteVector3dFieldFloat64 (
file, /*!< IN: file
handle */
"Hfield", /*!< IN: name of
dataset to write */
FieldstrengthHx, /*!< IN: X axis
data */
FieldstrengthHy, /*!< IN: Y axis
data */
FieldstrengthHz /*!< IN: Z axis
data */
);
H5Block3dSetFieldSpacing(file, "Hfield", stepHx, stepHy, stepHz);
H5Block3dSetFieldOrigin (file, "Hfield", sHx[0], sHy[0], sHz[0]);
}

H5WriteFileAttribFloat64 (
file, /*!< [in] Handle to
open file */
"Resonance Frequency(Hz)", /*!< [in] Name of
attribute */
&freq, /*!< [in] Array of
attribute values */
1 /*!< [in] Number of
array elements */
);
H5CloseFile(file);

std::cout << "Done bye ..." << std::endl;
std::cout << "--------------------------------------------------------"
<< std::endl;
}


Archive powered by MHonArc 2.6.19.

Top of Page