k-Wave

1. Anthony
   

   Member
   Joined: Apr '13
   Posts: 96

   Hi Brad,

   I would like to save only uz_non_staggered in my result file because don't need the other components and u_non_staggered is simply too big.

   I guess there is juste one or two lines to comment in the source code but I can't figure out precisely what to do. My hypothesis is that it should work if I comment the lines 705 to 714 in MatrixContainer.cpp but I am not familiar enough with the methods you use to write the HDF5 files...

   Best regards,
   Anthony

   PS: by the way, it was convenient to be able to save directly intensity in k-wave 1.0 (instead of pressure and velocity), I would love to see it back in 1.2 :-)

   Posted 9 years ago #

2. Jiri Jaros
   

   Developer
   Joined: Feb '12
   Posts: 115

   Hi Hypermoi,

   To only store uz_non_staggered you need to comment out a couple of lines at different places in the C code.

   This is not a clear solution

   1. MatrixContainer.cpp
   You have to inspect the `TMatrixContainer::AddMatricesIntoContainer() routine and remove these lines
   

   MatrixContainer[ux_shifted].SetAllValues(....);
   MatrixContainer[uy_shifted].SetAllValues(....);

   This removes the matrices from the Matrix container (prevents allocation of memory)

   2. TOutputStreamContainer::AddStreamsIntoContainer()
   Here, remove these two lines
   

   OutputStreamContainer[ux_shifted_sensor_raw] = CreateNewOutputStream(...)
   OutputStreamContainer[uy_shifted_sensor_raw] = CreateNewOutputStream(...)

   which will stop sampling the quantities.

   3. You need to disable calculating the x and y components

   Inspect TKSpaceFirstOrder3DSolver::Calculate_shifted_velocity() and remove the code for X and Y

   If it does not work, let me know and I'll send you a patch.

   I might also add a command line parameter to let you choose which components to sample, however we thought it would not be of any use :)

   Posted 9 years ago #

3. Jiri Jaros
   

   Developer
   Joined: Feb '12
   Posts: 115

   To answer the second question "Calculate intensity in the C code", I have to reply that it is impractical due to memory limitations.

   As you may know, p and u are staggered both in space and time to improve accuracy (Brad could explain that in detail). In k-Wave 1.0, we used a linear interpolation to get u at the same grid as p, however it had to low accuracy.

   In k-wave 1.1 we introduced fft-based shifts that improved accuracy. The down side is that to do so in space, you need all spatial data for ux, uy and uz even if you sample 1 single grid point! Thus to sample even a tiny region you need 6 more matrices! (speaking about a domain 512^3, we need additional 3GBs of RAM). Let me say tat This is the easy part.

   Now it starts to be tricky, to do the same shift in TIME we need the whole time history - over the sampling period. I talk here about tens to hundreds of GBs of RAM...:-)

   The only solution for now is to progressively store the non-staggered velocity to file, then terminate the simulation, load the data in a different order (1 gridpoint at time, but at all time steps), run the FFT shift and store the data back to file. - this is not implemented in the C code now, and you have to do that in Matlab (I agree that's not very convenient).

   Cheers
   Jiri

   Posted 9 years ago #

4. Anthony
   

   Member
   Joined: Apr '13
   Posts: 96

   Hi Jiri,

   Thanks a lot for your reply: the solution you describe seems to work according to my first tests.

   I understand the problem of memory and it is definitely a good reason not to allow for computing the intensity directly in the simulation code.

   You answer just raised a question in my mind: I thought that "u_unstaggered" was unstaggered both in space and time but your post let me think that it is only unstaggered in space, isnt'it?

   Actually, the reason why I only want to save uz is that I compute the integral of Iz over each XY plane to check the evolution of the power along the propagation (so I can check the increase in absorption due non linear effects around the focus for example...).

   Cheers,
   Anthony

   Posted 9 years ago #

5. Bradley Treeby
   

   Developer
   Joined: Mar '10
   Posts: 1,643

   Hi Anthony,

   You are correct - the output 'u_non_staggered' is only non-staggered in space, it is still staggered in time. To calculate intensity within the MATLAB code, this is then passed to the timeShift function before being multiplied by the pressure. You can take a look inside the private function kspaceFirstOrder_saveIntensity.m to see exactly what is done. Once we figure out an efficient and general way to do it, we'd like to put the calculate of intensity (and heat deposition) back inside the C++ code. I agree it would be more convenient this way.

   Brad.

   Posted 9 years ago #

6. Anthony
   

   Member
   Joined: Apr '13
   Posts: 96

   Hi Brad,

   Thanks for this important remark (and the pointer to the kspaceFirstOrder_saveIntensity.m function) :-) By the way, I hadn't seen it anywhere in the doc, maybe it would be worth adding it or mentioning it in the release notes, where you speak about "the additional input options".

   If you want to perform only heat deposition calculation, maybe you could only unstagger the last period as you basically only want to have a mean deposited power over one period. It would save a lot of memory...

   Anthony

   Posted 9 years ago #

7. Bradley Treeby
   

   Developer
   Joined: Mar '10
   Posts: 1,643

   Hi Anthony,

   Yes we did consider this for the calculation of I_avg or Q for CW fields. It would involve the user having to define the length of the acoustic period over which to average, and preferably defining the time discretisation to use an integer number of points per period. However, it would save a lot of memory as you say. Any other suggestions most welcome!

   Brad.

   Posted 9 years ago #

RSS feed for this topic

Reply

You must log in to post.