Setting up a first Clara2 simulation

[PrometheusPi]

This is a continuation of issue #93
Since the discussion changed from how to set up the environment on tianhe2 to how to set up Clara2 in general, I think a new issue should be used (to avid scrolling down through all the issues related to the module system on tinahe2)

This question was asked by @QJohn2017 in #93

---

In the beginning of settings.hpp you find the basic parameters:

• omega_max is the maximum frequency of radiation you want to resolve
  so your frequency range will start at 0.0 1/s and go to omega_max 1/s

• theta_max is the maximum opening angle of your detector in degree
  Setting the detector is not yet generic and requires some code adjustments we will discuss later

• N_spectrum is the number of frequencies you want to sample between 0.0 1/s and go to omega_max 1/s

• N_theta is the number of angle values between 0 degree and theta_max degree you want to sample

• N_phi is the number of angles you want to sample orthogonally to theta

• N_trace is the number of trajectories provided that you want to calculate the spectrum for.

• fft_length_factor is a factor to increase the sampling of the trajectories via linear interpolation. This is essential when you undergo highly-nonlinear Thomson scattering for only a brief period of time compared to the entire duration of the trajectory (see https://doi.org/10.5281/zenodo.843510 section 4.1.3)

• ascii_output if true returns the spectra as text file, if false returns the spectra as binary file

• N_char_filename the number of characters that are needed to give the trajectory path. 256 are usually sufficient, but if you have long directory names, more might be needed

• traceFileTemplate the template for locating your input trajectories with [C-style])http://www.cplusplus.com/reference/cstdio/printf/) replacements

• outputFileTemplate the template for the output spectra files with C-style replacements

[PrometheusPi]

In the input section of settings.hpp the data in the trajectory file is defined. By default a GPT output trajectory is assumed that looks like:

x_position[meter]    y_position[meter]    z_position[meter]   beta_x[unitless]  beta_y[unitless]   beta_z[unitless]   time[seconds] 

with beta being the velocity normalized to the speed of light and beta_{x,y,z} being its vector components.

Each column contains the value specified above.
Each row contains a new point of the trajectory with an increased time.

This results in the following input settings:

namespace input
{
  const unsigned int index_time = 6; /* column id of time */

  const unsigned int index_pos_x = 0; /* column id of x */
  const unsigned int index_pos_y = 1; /* column id of y */
  const unsigned int index_pos_z = 2; /* column id of z */

  const unsigned int index_beta_x = 3; /* column id of beta_x */
  const unsigned int index_beta_y = 4; /* column id of beta_y */
  const unsigned int index_beta_z = 5; /* column id of beta_z */

  /* unit conversion to SI units */
  const double convert_time = 1.0; /* s -> s */
  const double convert_pos = 1.0; /* m -> m */
  const double convert_beta = 1.0; /* None -> None */
}

If you for example prefere cgs units and the trajectory you have is in units of

• time: seconds (s)
• position: centimeter (cm)
• speed: cm/s
  with the above order, the input settings would looks like

namespace input
{
  const unsigned int index_time = 0; /* column id of time */

  const unsigned int index_pos_x = 1; /* column id of x */
  const unsigned int index_pos_y = 2; /* column id of y */
  const unsigned int index_pos_z = 3; /* column id of z */

  const unsigned int index_beta_x = 4; /* column id of beta_x */
  const unsigned int index_beta_y = 5; /* column id of beta_y */
  const unsigned int index_beta_z = 6; /* column id of beta_z */

  /* unit conversion to SI units */
  const double convert_time = 1.0; /* s -> s */
  const double convert_pos = 0.01; /* cm -> m */
  const double convert_beta = 1.0/2.9978e10; /* cm/s -> None */
}

[PrometheusPi]

If you have a more complex input set, for example if you use momentum values instead of velocities, you should have a look at run_trough_data.hpp. There you could add more complex (non-linear) conversions for your specific trajectory data format.

[PrometheusPi]

The last part of settings.hpp is for the post- processing via process defined here. It collects the radiation of all trajectories and adds them up, to generate a total (incoherent) spectra.

Nothing needs to be changed there.

[TheresaBruemmer]

This last part is what does not work in my case. The single trace spectra are not added up. But I would expect a "could not write output file" ... but I get no error at all...

[PrometheusPi]

Regarding the observation points:

The number of total directions is the product of N_theta and N_phi defined in settings.hpp.

Their individual values in degree are defined here. Be default phi has the values 0 degree and 90 degree.

These values are handed to the single_direction.cpp code for each combination individually here.

In the single_direction.cpp code, the actual unit vector pointing towards the observer is defined here.
First of all, all angles are converted to radiant. And then the vector is computed.
If you want to change the observation direction, please adjust the code of the vector. By default, the vector is

looking_vector = R_vec(std::cos(angle_theta) ,
                         std::sin(angle_theta) * std::cos(angle_phi)  ,
std::sin(angle_theta) * std::sin(angle_phi) );

so for phi = 0 degree this is (cos(theta), sin(theta), 0) while for phi = 90 degreethis is (cos(theta), 0, sin(theta)).
This observer points for theta=0 toward the +x direction (1,0,0) regardless of phi.

[PrometheusPi]

@TheresaBruemmer In order to perform the post-processing you need to run

./process_data

after the parallel job has finished.

The parallel job should have created 3 types of files:

• my_output.txt... for output information on each trajectory
• my_error.txt... for error information on each trajectory
• my_spectrum_trace......dat for the individual spectra for each trajectory

Did you run ./process_data?

Only this last step creates the my_spectrum_all***.txt files.

[TheresaBruemmer]

Nice, thank you. It works!

[PrometheusPi]

In a future version, where no array jobs are supported any more, this step can be executed by clara2 directly. But with array job support, we need a separate executable for this step, because the individual array jobs are not aware of each other. With MPI only, this would not be required.

[PrometheusPi]

@QJohn2017 Did setting up your simulation work for you?

[QJohn2017]

@PrometheusPi Sorry that I'm busy these days and have no time to set my simulation, i will have a try to create the trajectory files and set the simulation again tonight. My previous trajectory files are generated by the Vsim code and only has the position information of electrons that are being tracked.

[QJohn2017]

@PrometheusPi , I have set up the first simulation successfully, when I input the command./process_data, there are two output files named my_spectrum_all_000.dat and my_spectrum_all_001.dat created. But I don't know how these two files stand for what and how to analyse them

[PrometheusPi]

@QJohn2017 No problem - thanks for the reply.
Great to hear that running Clara2 worked with your data 👍. That are wonderful news.

The two files ./process_data creates are the directionally resolved spectra for both phi values you used. So in case you were using the default setup my_spectrum_all_000.dat is the phi=0° spectrum and my_spectrum_all_001.dat is the phi=90° spectrum.
Both files are text files and are a matrix with columns describing each frequency omega values and the rows describing your theta angle values.

You can load these files directly into your favorite data anaysis/visualization framework like Mathematica, Python (numpy, matplotlib, ...), MatLab, etc.

There is also a (python based) tool provided in clara2/tools/plotRadiation.
If you run

./clara2/tools/plotRadiation my_spectrum_all_000.dat

it will plot the directionally resolved spectrum for phi=0° to your screen (if logged in with X-Window via ssh -X ...). Without arguments this uses default values for theta and omega. There is an extensive help that allows you to specify these values as well. You can access this help by running ./plotRadiation -h. If you for example computed radiation from 0 to 100 Hz, and your theta angle goes from -10° to +30°, and you want to produce a pdf instead of a window on your screen, you could run:

./clara2/tools/plotRadiation --dataExtend 0 100 -10 30  --labelOmega "$\omega\ [Hz]$" --pdf my_spectrum_all_000.dat 

If you have questions how to load these output files into python, please feel free to ask me. With all other post processing tools I am unfortunately not that familiar.

@QJohn2017 Did the spectra came out as you expected?