NanoParticles distribution

Dear Geant4 users, I face a problem in simulation of Nanoparticles in my code, would any one can give me a clue how I can tackle this problem? the position of the Nanoparticles in the Geant4 code is as shown in snapshot below.

The above snapshot code is working for NPs = 107 but when I increased the number to 1010 and above it couldn’t simulate, would any one helped me how I can manage this problem please?

  1. Are you sure “NP” (the number of nanoparticles) i supposed to be a double, rather than an integer?

  2. You don’t show any code where that variable is being used, which makes it hard for anyone other than you to know how it’s supposed to be used.

  3. Please don’t use screenshots for code. Copy and paste the code as text, using a “block” delimiter of three backquotes at the beginning and end:

    void myFunction(G4int someArg) {
      G4cout << someArg << " Hello, world!" << G4endl;

which will be displayed as:

void myFunction(G4int someArg) {
  G4cout << someArg << " Hello, world!" << G4endl;

Thank you very much for response and Sorry, for the snapshot post. Yes you are true for NPs G4int or G4intLong is the true representation, how ever, the number of of the NPs are beyond it. When I run it gives error, but I changed to double it works and generate the nano particles. My challenge is when i increase the number beyond 10*7 it takes weeks and even doesn’t show any result. So is there any mechanism that makes my simulation fast? here is the code. To distribute the NPs I used std::vector and G4UniformRand
// Cylindrical Water Phantom

G4double pRmin =;
G4double pRmax = 3.25
G4double pHeight = 7.5cm;
G4double pSPhi = 0.deg;
G4double pDPhi = 360.deg;
solidTube = new G4Tubs(“solidTube”, pRmin, pRmax, pHeight, pSPhi, pDPhi);
logicTube = new G4LogicalVolume(solidTube, water,“logicTube”);
physTube = new G4PVPlacement(0, G4ThreeVector(0, 0, 0),
logicTube, “physTube”, logicPhantom, false, 0, checkOverlaps);

G4VisAttributes* solidTubeVisAtt = new G4VisAttributes(brown);
solidTubeVisAtt → SetVisibility(true);
solidTubeVisAtt → SetForceWireframe(true);
logicTube → SetVisAttributes(solidTubeVisAtt);

// generating positions inside the cylindrical phantom

G4double NP = 1.8e10;       // Number of Sherical gold Nanoparticles in (18mg/g)
G4double r =  5e-6*cm;       //radius of the Sherical Nanoparticles
G4double DMIN = 2*r;         // Minimum distance between two Sherical Nanoparticles

G4ThreeVector fmin, fmax;
solidTube->BoundingLimits(fmin,fmax); // get bounding box of the cylinder
std::vector position;
while(position.size() < NP)

if you place more than G4INT_MAX (~2e9) volumes, the instanceID of type G4int for the physical volume will overflow. don’t know if that necessarily causes problems, but apparently it very likely does so :slight_smile:

also consider the probability to generate a free spot for the next NP to place. if the volume is already critically filled, it will take significantly longer and longer to randomly generate a non-overlapping location.

Such a huge number of objects as 1.8*10^10 nanoparticles cannot be generated in one volume. Just filling a std::vector with nanoparticle positions would require half a terabyte of operating memory and hundreds of years of CPU time:

  • One position (G4ThreeVector) takes 24 bytes, 1.8*10^10 positions will take 432 Gb

  • Current algorithm to generate NP positions performs more than NP^2/2 comparisons. Assuming that NP/2 = 0.9*10^10 comparisons can be done in 1 sec (actually more), it means that all comparisons will take 1.8*10^10 secs = 570 years

So, you should use much smaller volumes. For example, reduce the size of solidTube:

pRmax = 3.25*cm / 25; // 0.13*cm
pHeight = 7.5*cm / 25; // 0.3*cm
NP = 1.152e6;

Alternatively, you can try to construct the original volume as a set of sub-volumes (multiple replication of a very small volume).

Surprising, Thank you very much evc for your clarification. I will try both of your ways. In the 1st case when you reduce the pRmax and pHeight of the cylinder by 25, you put the NP = 1.152e6; Would you please explain how you determine the NP value = 1.152e6?. I mean from where comes from this 1.152e6 number?

Both, the radius and height of the cylinder are reduced in 25 times, it means that the volume of the cylinder is reduced in 25*25*25 = 15625 times. To keep the density of nanoparticles in the volume unchanged, the number of nanoparticles should also be reduced in 15625 times:

1.8*10^10 / (25*25*25) = 1.152*10^6

Thank you evc once again, as you said it is better to reduce the size of the cylinder. “Alternatively, you can try to construct the original volume as a set of sub-volumes (multiple replication of a very small volume)” I need also to apply this idea in the code below.**

G4Box* tumor = new G4Box(“tumor”, 0.5cm, 0.5cm, 0.5cm);
logicTumor = new G4LogicalVolume(tumor, softTissue, “logicTumor”);
G4VPhysicalVolume* physTumor = new G4PVPlacement(0, G4ThreeVector(0, 1*cm, 0), logicTumor, “physTumor”, logicTube, false, 0, checkOverlaps);

Here, I need to divide the tumor volume of 1cc in 125000 parts, that means in each part we will have 1.44*10^4 NPs. i.e each side of the cube is divided by 50.

The minimum and maximum value of the cube are[-0.5, 0.5) = G4ThreeVector fmin, fmax;

To distribute the NPs along the x-axis without dividing the cube we can use

G4double x = fmin.x() + (fmax.x() -fmin.x())*G4UniformRand();

but after dividing the into 125000 parts how we can distribute the NPs?

There are three logical volumes: big box, small box, and nanoparticle.

You fill the small box with randomly distributed 1.44*10^4 nanoparticles using an algorithm very similar to the one that was used for the cylinder.

Then, in three nested for loops, you fill the big box with 50x50x50 small boxes using corresponding translation vectors (G4ThreeVector).

Dear evc really appreciate for all you thoughtful suggestions. Based on the quoted above recommendation I wrote the three nested loops as below but I couldnt succeeded. Would you please suggest me what mistake I made?
//// Voxel
G4Box* voxel = new G4Box(“voxel”, 0.1mm, 0.1mm, 0.1mm);
logicVoxel = new G4LogicalVolume(voxel, softTissue, “logicVoxel”);
G4int N =50;
for(G4int j = 0; j < N; j++)


  for(G4int m = 0; m < N-2; m++)
  for(G4int n = 0; n < N; n++)

  G4VPhysicalVolume* physVoxel = new G4PVPlacement(0, G4ThreeVector(-0.5*cm+(j+0.5*cm)/N), 1*cm+(m+2)/N),-0.5*cm             +(n+0.5*cm)/N)), 0), logicVoxel, "physVoxel", logicTumor, false, 0, checkOverlaps);

*The minimum and maximum value of the cube are[-0.5, 0.5)

Try the code below in visual mode. If the image looks Ok, then you can try to replace the parameters in the code with the ones you need.

  // Main parameters                                                                                                    
  G4double sizeTumor = 10*mm;
  G4double R = 0.05*mm; // radius of particle                                                                           
  G4int Nv = 7;         // number of voxels = Nv * Nv * Nv                                                              
  G4int Np = 100;       // number of particles in voxel                                                                 
  G4double sizeVoxel = sizeTumor / Nv;

  // Materials                                                                                                          
  G4Material* tumor_mat = nist->FindOrBuildMaterial("G4_WATER");
  G4Material* voxel_mat = tumor_mat;
  G4Material* nano_mat = nist->FindOrBuildMaterial("G4_Au");

  // Make logical volumes                                                                                               
  auto solidTumor = new G4Box("Tumor", 0.5*sizeTumor, 0.5*sizeTumor, 0.5*sizeTumor);
  auto logicTumor = new G4LogicalVolume(solidTumor, tumor_mat, "Tumor_LV");

  auto solidVoxel = new G4Box("Voxel", 0.5*sizeVoxel, 0.5*sizeVoxel, 0.5*sizeVoxel);
  auto logicVoxel = new G4LogicalVolume(solidVoxel, voxel_mat, "Voxel_LV");

  auto solidNano = new G4Orb("Nano", R);
  auto logicNano = new G4LogicalVolume(solidNano, nano_mat, "Nano_LV");

  // Fill Voxel with nanoparticles                                                                                      
  std::vector<G4ThreeVector> pos;
  G4double DMIN = 2*R;
  while (pos.size() < Np) {
    // generate position                                                                                                
    G4double x = (0.5*sizeVoxel - R)*(2.*G4UniformRand() - 1.);
    G4double y = (0.5*sizeVoxel - R)*(2.*G4UniformRand() - 1.);
    G4double z = (0.5*sizeVoxel - R)*(2.*G4UniformRand() - 1.);
    G4ThreeVector p(x, y, z);
    // check position                                                                                                   
    G4int size = pos.size();
    G4int k = 0;
    for ( ; k < size; ++k) { if ((pos[k] - p).mag() < DMIN) break; }
    // add position to the list                                                                                         
    if (k == size) pos.push_back(p);
  // Place nanoparticles at generated positions                                                                         
  for (G4int i = 0; i < Np; ++i) {
    new G4PVPlacement(nullptr,pos[i],logicNano, "Nano_PV", logicVoxel, false,0,false);

  // Fill Tumor with Voxels  
  G4double dd = 0.5*sizeVoxel;
  for (G4int ix = 0; ix < Nv; ++ix) {
    G4double x = -(0.5*sizeTumor - dd) + 2*ix*dd;
    for (G4int iy = 0; iy < Nv; ++iy) {
      G4double y = -(0.5*sizeTumor - dd) + 2*iy*dd;
      for (G4int iz = 0; iz < Nv; ++iz) {
        G4double z = -(0.5*sizeTumor - dd) + 2*iz*dd;
        new G4PVPlacement(nullptr,G4ThreeVector(x,y,z),logicVoxel,"Voxel_PV",logicTumor, false,0,false);

This topic was automatically closed 7 days after the last reply. New replies are no longer allowed.