ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4QAOLowEnergyLoss.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4QAOLowEnergyLoss.cc
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 //
27 // -------------------------------------------------------------
28 // GEANT 4 class implementation file
29 //
30 // History: New Implementation
31 //
32 // ---------- G4QAOLowEnergyLoss physics process -------
33 // by Stephane Chauvie, 5 May 2000
34 // Modified:
35 //
36 // 24/05/2000 MGP Modified to remove compilation warnings on Linux and DEC
37 // Introduced sizes of L0, L1, L2 arrays
38 // 23/05/2000 MGP Made compliant to design
39 // 02/08/2000 V.Ivanchenko Clean up according new design
40 // 16/09/2000 S. Chauvie Oscillator for all materials
41 // 03/10/2000 V.Ivanchenko CodeWizard clean up
42 // 05/11/2000 V.Ivanchenko "Aluminum" - correct name, end of cycle
43 // over shells, and two bugs from previous edition
44 // 10/05/2001 V.Ivanchenko Clean up againist Linux compilation with -Wall
45 // 13/05/2001 S. Chauvie corrected bugs
46 // 01/06/2001 V.Ivanchenko replace names by Z, change the validity range
47 // from 50 keV to 5 KeV and change sign of the
48 // Barkas term
49 // 4/06/2001 S. Chauvie Corrected small bugs
50 //
51 // ************************************************************
52 // It is the Quantal Harmonic Oscillator Model for energy loss
53 // of slow antiproton
54 // ************************************************************
55 // --------------------------------------------------------------
56 
57 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
58 
59 #include "G4QAOLowEnergyLoss.hh"
60 #include "G4PhysicalConstants.hh"
61 #include "G4SystemOfUnits.hh"
62 #include "G4DynamicParticle.hh"
63 #include "G4Material.hh"
64 #include "G4ParticleDefinition.hh"
65 #include "G4AntiProton.hh"
66 #include "G4Exp.hh"
67 
68 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
69 
71  : G4VLowEnergyModel(name)
72 {
74  sizeL0 = 67;
75  sizeL1 = 22;
76  sizeL2 = 14;
77 }
78 
79 
81 {;}
82 
83 
85  const G4Material* ) const
86 {
87  return 2.0*MeV ;
88 }
89 
90 
92  const G4Material* ) const
93 {
94  // return 50.0*keV ;
95  return 5.0*keV ;
96 }
97 
98 
100 {
101  return 2.0*MeV ;
102 }
103 
105 {
106  // return 50.0*keV ;
107  return 5.0*keV ;
108 }
109 
110 
112  const G4Material* material) const
113 {
114  G4bool isInCharge = false;
115 
116  G4bool hasMaterial = false;
117 
118  if (material->GetNumberOfElements() == 1) hasMaterial = true;
119 
120  if ((particle->GetDefinition()) == (G4AntiProton::AntiProtonDefinition())
121  && hasMaterial) isInCharge = true;
122 
123  return isInCharge;
124 
125 }
126 
128  const G4Material* material) const
129 {
130 
131  G4bool isInCharge = false;
132 
133  G4bool hasMaterial = false;
134 
135  if (material->GetNumberOfElements() == 1) hasMaterial = true;
136 
137 
138  if (aParticle == (G4AntiProton::AntiProtonDefinition())
139  && hasMaterial) isInCharge = true;
140 
141  return isInCharge;
142 
143 }
144 
145 
147  const G4Material* material)
148 {
149  G4double zParticle = (G4int)(particle->GetCharge())/eplus;
150 
151  G4double energy = particle->GetKineticEnergy() ;
152  G4double eloss = EnergyLoss(material,energy,zParticle) ;
153 
154  return eloss ;
155 }
156 
157 
159  const G4Material* material,
160  G4double kineticEnergy)
161 {
162  G4double zParticle = (aParticle->GetPDGCharge())/eplus;
163 
164  G4double eloss = EnergyLoss(material,kineticEnergy,zParticle) ;
165 
166  return eloss ;
167 }
168 
169 
170 
172  G4double kineticEnergy,
173  G4double zParticle) const
174 {
175  G4int nbOfShell = GetNumberOfShell(material);
176  if(nbOfShell < 1) nbOfShell = 1;
177 
178  //G4int iz = (G4int)(material->GetZ());
179  //if(iz < 1) iz = 1;
180  //else if(iz > 100) iz = 100;
181  //G4int nbOfShell = fNumberOfShells[iz];
182 
183  /*
184  if(material->GetName()=="Graphite"){
185  G4cout << " E(MeV)= " << kineticEnergy/MeV << " n= "
186  << nbOfShell
187  << " for " << material->GetZ() << G4endl;
188  }
189  */
190  G4double dedx=0;
191 
192  G4double v = c_light * std::sqrt( 2.0 * kineticEnergy / proton_mass_c2 );
193  G4double coeff = twopi * proton_mass_c2 *
194  (material-> GetTotNbOfElectPerVolume()) /
196  G4double fBetheVelocity = fine_structure_const * c_light / v;
198  kineticEnergy ;
199 
200  //G4double beta = std::sqrt( 2.0 * kineticEnergy / proton_mass_c2 );
201  //G4double fBetheVelocity = std::sqrt( 2.0 * 25.0 * keV / proton_mass_c2 )/beta;
202  //G4double coeff= twopi_mc2_rcl2*(material->GetElectronDensity())/(beta*beta);
203 
204 
205  G4double l0Term = 0, l1Term = 0, l2Term = 0;
206 
207  for (G4int nos = 0 ; nos < nbOfShell ; nos++){
208 
209  G4double l0 = 0, l1 = 0, l2 = 0;
210  G4double NormalizedEnergy = ( 2.0 * electron_mass_c2 * v * v ) /
211  ( c_squared * GetShellEnergy(material,nos) );
212  //G4double NormalizedEnergy = ( 2.0 * electron_mass_c2 * beta * beta ) /
213  // GetShellEnergy(material,nos);
214  G4double shStrength = GetShellStrength(material,nos);
215 
216  l0 = GetL0(NormalizedEnergy);
217  l0Term += shStrength * l0;
218 
219  l1 = GetL1(NormalizedEnergy);
220  l1Term += shStrength * l1;
221 
222  l2 = GetL2(NormalizedEnergy);
223  l2Term += shStrength * l2;
224 
225  /*
226  if(material->GetName()=="Graphite"){
227  G4cout << nos << ". "
228  << " E(MeV)= " << kineticEnergy/MeV
229  << " normE= " << NormalizedEnergy
230  << " sh en= " << GetShellEnergy(material,nos)
231  << " str= " << shStrength
232  << " v0/v= " << fBetheVelocity
233  << " l0= " << l0Term
234  << " l1= " << l1Term
235  << " l2= " << l2Term
236  << G4endl;
237  }
238  */
239  }
240 
241  dedx = coeff * zParticle * zParticle * (l0Term
242  + zParticle * fBetheVelocity * l1Term
243  + zParticle * zParticle * fBetheVelocity * fBetheVelocity * l2Term);
244 
245  // G4cout << " E(MeV)= " << kineticEnergy/MeV
246  // << " dedx(Mev/mm)= " << dedx*mm/MeV << G4endl;
247 
248  return dedx ;
249 
250 }
251 
252 
254 {
255  // Set return value from table
256  G4int Z = (G4int)(material->GetZ());
257  G4int nShell = 0;
258 
259  // Set return value if in material available from Aahrus
260  for(G4int i=0; i<numberOfMaterials; i++) {
261 
262  if(materialAvailable[i] == Z){
263  nShell = nbofShellForMaterial[i];
264  break;
265  }
266  else nShell = fNumberOfShells[Z];
267  }
268 
269  return nShell;
270 }
271 
272 
273 
275  G4int nbOfTheShell) const
276 {
277  //
278  G4double shellEnergy = alShellEnergy[0];
279 
280  if(material->GetZ() == 13) shellEnergy = alShellEnergy[nbOfTheShell];
281  else if(material->GetZ() == 14)shellEnergy = siShellEnergy[nbOfTheShell];
282  else if(material->GetZ() == 29)shellEnergy = cuShellEnergy[nbOfTheShell];
283  else if(material->GetZ() == 73)shellEnergy = taShellEnergy[nbOfTheShell];
284  else if(material->GetZ() == 79)shellEnergy = auShellEnergy[nbOfTheShell];
285  else if(material->GetZ() == 78)shellEnergy = ptShellEnergy[nbOfTheShell];
286  else if (material->GetNumberOfElements() == 1)
287  shellEnergy = GetOscillatorEnergy(material, nbOfTheShell);
288  else G4cout << "WARNING - G4QAOLowEnergyLoss::GetShellEnergy - "
289  << "The model is not available for "
290  << material->GetName()
291  << G4endl;
292 
293  return shellEnergy;
294 }
295 
296 
298  G4int nbOfTheShell) const
299 {
300 
301  const G4Element* element = material->GetElement(0);
302 
303  G4int Z = (G4int)(element->GetZ());
304 
305 G4double squaredPlasmonEnergy = 28.816 * 28.816 * 1e-6
306  * material->GetDensity()/g/cm3
307  * (Z/element->GetN()) ;
308 //G4double squaredPlasmonEnergy = 28.16 * 28.16 * 1e-6
309 // * (material->GetDensity()) * (cm3/g)
310 // * Z / (element->GetN()) ;
311 
312  G4double plasmonTerm = 0.66667 * GetOccupationNumber(Z,nbOfTheShell)
313  * squaredPlasmonEnergy / (Z*Z) ;
314 
315  G4double ionTerm = G4Exp(0.5) * (element->GetAtomicShell(nbOfTheShell)) ;
316 
317  ionTerm = ionTerm*ionTerm ;
318 
319  G4double oscShellEnergy = std::sqrt( ionTerm + plasmonTerm );
320 
321 /* if(material->GetName()=="Graphite"){
322  G4cout << "\t" << Z
323  << "\t" << nbOfTheShell
324  << "\t" << squaredPlasmonEnergy
325  << "\t" << plasmonTerm
326  << "\t" << ionTerm
327  << "\t" << GetOccupationNumber(Z,nbOfTheShell)
328  << "\t" << element->GetAtomicShell(nbOfTheShell)
329  << "\t" << oscShellEnergy << G4endl;}
330 */
331  return oscShellEnergy;
332 }
333 
335  G4int nbOfTheShell) const
336 {
337  G4double shellStrength = alShellStrength[0];
338 
339  if(material->GetZ() == 13) shellStrength = alShellStrength[nbOfTheShell];
340  else if(material->GetZ() == 14)shellStrength =siShellStrength[nbOfTheShell];
341  else if(material->GetZ() == 29)shellStrength =cuShellStrength[nbOfTheShell];
342  else if(material->GetZ() == 73)shellStrength =taShellStrength[nbOfTheShell];
343  else if(material->GetZ() == 79)shellStrength =auShellStrength[nbOfTheShell];
344  else if(material->GetZ() == 78)shellStrength =ptShellStrength[nbOfTheShell];
345  else if (material->GetNumberOfElements() == 1){
346  G4int Z = (G4int)(material->GetZ());
347  shellStrength = GetOccupationNumber(Z,nbOfTheShell) / Z ;}
348  else G4cout << "WARNING - G4QAOLowEnergyLoss::GetShellEnergy - "
349  << "The model is not available for "
350  << material->GetName()
351  << G4endl;
352 
353  return shellStrength;
354 }
355 
357 {
358 
359  G4int indice = ShellNb ;
360  for (G4int z = 1 ; z < Z ; z++) {indice += fNumberOfShells[z];}
361 
362  return nbOfElectronPerSubShell[indice+1];
363 }
364 
365 
367 {
368  G4int n;
369 
370  for(n = 0; n < sizeL0; n++) {
371  if( normEnergy < L0[n][0] ) break;
372  }
373  if(0 == n) n = 1 ;
374  if(n >= sizeL0) n = sizeL0 - 1 ;
375 
376  G4double l0 = L0[n][1];
377  G4double l0p = L0[n-1][1];
378  G4double bethe = l0p + (l0 - l0p) * ( normEnergy - L0[n-1][0]) /
379  (L0[n][0] - L0[n-1][0]);
380  return bethe ;
381 
382 }
383 
385 {
386  G4int n;
387 
388  for(n = 0; n < sizeL1; n++) {
389  if( normEnergy < L1[n][0] ) break;
390  }
391  if(0 == n) n = 1 ;
392  if(n >= sizeL1) n = sizeL1 - 1 ;
393 
394  G4double l1 = L1[n][1];
395  G4double l1p = L1[n-1][1];
396  G4double barkas= l1p + (l1 - l1p) * ( normEnergy - L1[n-1][0]) /
397  (L1[n][0] - L1[n-1][0]);
398 
399  return barkas;
400 
401 }
402 
403 
405 {
406  G4int n;
407  for(n = 0; n < sizeL2; n++) {
408  if( normEnergy < L2[n][0] ) break;
409  }
410  if(0 == n) n = 1 ;
411  if(n >= sizeL2) n = sizeL2 - 1 ;
412 
413  G4double l2 = L2[n][1];
414  G4double l2p = L2[n-1][1];
415  G4double bloch = l2p + (l2 - l2p) * ( normEnergy - L2[n-1][0]) /
416  (L2[n][0] - L2[n-1][0]);
417 
418  return bloch;
419 }
420 
421 
422 const G4int G4QAOLowEnergyLoss::materialAvailable[6] = {13,14,29,73,79,78};
423  /*
424  "Aluminum",
425  "Silicon",
426  "Copper",
427  "Tantalum",
428  "Gold",
429  "Platinum"};
430  */
431 const G4int G4QAOLowEnergyLoss::nbofShellForMaterial[6] = {3,3,4,6,6,6 };
432 
433 const G4double G4QAOLowEnergyLoss::alShellEnergy[3] ={ 2795e-6, 202e-6, 16.9e-6};
434 const G4double G4QAOLowEnergyLoss::alShellStrength[3]={ 0.1349, 0.6387, 0.2264};
435 const G4double G4QAOLowEnergyLoss::siShellEnergy[3] ={ 3179e-6, 249e-6, 20.3e-6 };
436 const G4double G4QAOLowEnergyLoss::siShellStrength[3]={ 0.1222, 0.5972, 0.2806};
437 const G4double G4QAOLowEnergyLoss::cuShellEnergy[4] ={ 16931e-6, 1930e-6, 199e-6, 39.6e-6};
438 const G4double G4QAOLowEnergyLoss::cuShellStrength[4]={ 0.0505, 0.2561, 0.4913, 0.2021};
439 const G4double G4QAOLowEnergyLoss::taShellEnergy[6] ={ 88926e-6, 18012e-6, 3210e-6, 575e-6, 108.7e-6, 30.8e-6};
440 const G4double G4QAOLowEnergyLoss::taShellStrength[6]={ 0.0126, 0.0896, 0.2599, 0.3413, 0.2057, 0.0908};
441 const G4double G4QAOLowEnergyLoss::auShellEnergy[6]={ 96235e-6, 25918e-6, 4116e-6, 599e-6, 87.3e-6, 36.9e-6};
442 const G4double G4QAOLowEnergyLoss::auShellStrength[6]={ 0.0139, 0.0803, 0.2473, 0.423, 0.1124, 0.1231};
443 const G4double G4QAOLowEnergyLoss::ptShellEnergy[6]={ 95017e-6, 25590e-6, 4063e-6, 576e-6, 81.9e-6, 31.4e-6};
444 const G4double G4QAOLowEnergyLoss::ptShellStrength[6]={ 0.0129, 0.0745, 0.2295, 0.4627, 0.1324, 0.0879};
445 
446 
447 const G4double G4QAOLowEnergyLoss::L0[67][2] =
448 {
449  {0.00, 0.000001},
450  {0.10, 0.000001},
451  {0.12, 0.00001},
452  {0.14, 0.00005},
453  {0.16, 0.00014},
454  {0.18, 0.00030},
455  {0.20, 0.00057},
456  {0.25, 0.00189},
457  {0.30, 0.00429},
458  {0.35, 0.00784},
459  {0.40, 0.01248},
460  {0.45, 0.01811},
461  {0.50, 0.02462},
462  {0.60, 0.03980},
463  {0.70, 0.05731},
464  {0.80, 0.07662},
465  {0.90, 0.09733},
466  {1.00, 0.11916},
467  {1.20, 0.16532},
468  {1.40, 0.21376},
469  {1.60, 0.26362},
470  {1.80, 0.31428},
471  {2.00, 0.36532},
472  {2.50, 0.49272},
473  {3.00, 0.61765},
474  {3.50, 0.73863},
475  {4.00, 0.85496},
476  {4.50, 0.96634},
477  {5.00, 1.07272},
478  {6.00, 1.27086},
479  {7.00, 1.45075},
480  {8.00, 1.61412},
481  {9.00, 1.76277},
482  {10.00, 1.89836},
483  {12.00, 2.13625},
484  {14.00, 2.33787},
485  {16.00, 2.51093},
486  {18.00, 2.66134},
487  {20.00, 2.79358},
488  {25.00, 3.06539},
489  {30.00, 3.27902},
490  {35.00, 3.45430},
491  {40.00, 3.60281},
492  {45.00, 3.73167},
493  {50.00, 3.84555},
494  {60.00, 4.04011},
495  {70.00, 4.20264},
496  {80.00, 4.34229},
497  {90.00, 4.46474},
498  {100.00, 4.57378},
499  {120.00, 4.76155},
500  {140.00, 4.91953},
501  {160.00, 5.05590},
502  {180.00, 5.17588},
503  {200.00, 5.28299},
504  {250.00, 5.50925},
505  {300.00, 5.69364},
506  {350.00, 5.84926},
507  {400.00, 5.98388},
508  {450.00, 6.10252},
509  {500.00, 6.20856},
510  {600.00, 6.39189},
511  {700.00, 6.54677},
512  {800.00, 6.68084},
513  {900.00, 6.79905},
514  {1000.00, 6.90474}
515 };
516 
517 const G4double G4QAOLowEnergyLoss::L1[22][2] =
518 {
519  {0.00, -0.000001},
520  {0.10, -0.00001},
521  {0.20, -0.00049},
522  {0.30, -0.00084},
523  {0.40, 0.00085},
524  {0.50, 0.00519},
525  {0.60, 0.01198},
526  {0.70, 0.02074},
527  {0.80, 0.03133},
528  {0.90, 0.04369},
529  {1.00, 0.06035},
530  {2.00, 0.24023},
531  {3.00, 0.44284},
532  {4.00, 0.62012},
533  {5.00, 0.77031},
534  {6.00, 0.90390},
535  {7.00, 1.02705},
536  {8.00, 1.10867},
537  {9.00, 1.17546},
538  {10.00, 1.21599},
539  {15.00, 1.24349},
540  {20.00, 1.16752}
541 };
542 
543 const G4double G4QAOLowEnergyLoss::L2[14][2] =
544 {
545  {0.00, 0.000001},
546  {0.10, 0.00001},
547  {0.20, 0.00000},
548  {0.40, -0.00120},
549  {0.60, -0.00036},
550  {0.80, 0.00372},
551  {1.00, 0.01298},
552  {2.00, 0.08296},
553  {4.00, 0.21953},
554  {6.00, 0.23903},
555  {8.00, 0.20893},
556  {10.00, 0.10879},
557  {20.00, -0.88409},
558  {40.00, -1.13902}
559 };
560 
562 {
563  0, // consistency with G4AtomicShells
564  1,//------ H
565  2,//------ He
566  2, 1,//------ Li
567  2, 2,//------ Be
568  2, 2, 1,//------ B
569  2, 2, 2,//------ C
570  2, 2, 2, 1,//------ N
571  2, 2, 2, 2,//------ O
572  2, 2, 5,//------ F
573  2, 2, 2, 4,//------ Ne
574  2, 2, 2, 4, 1,//------ Na
575  2, 2, 2, 4, 2,//------ Mg
576  2, 2, 2, 4, 2, 1,//------ Al
577  2, 2, 2, 4, 2, 2,//------ Si
578  2, 2, 2, 4, 2, 3,//------ P
579  2, 2, 2, 4, 2, 4,//------
580  2, 2, 2, 4, 2, 5,//------
581  2, 2, 2, 4, 2, 2, 4,//------
582  2, 2, 2, 4, 2, 2, 4, 1,//------
583  2, 2, 2, 4, 2, 2, 4, 2,//------
584  2, 2, 2, 4, 2, 2, 4, 1, 2,//------
585  2, 2, 2, 4, 2, 2, 4, 2, 2,//------
586  2, 2, 2, 4, 2, 2, 4, 3, 2,//------
587  2, 2, 2, 4, 2, 2, 4, 4, 2,//------
588  2, 2, 2, 4, 2, 2, 4, 5, 2,//------
589  2, 2, 2, 4, 2, 2, 4, 6, 2,//------
590  2, 2, 2, 4, 2, 2, 4, 7, 2,//------
591  2, 2, 2, 4, 2, 2, 4, 4, 4, 2,//------
592  2, 2, 2, 4, 2, 2, 4, 4, 5, 2,//------
593  2, 2, 2, 4, 2, 2, 4, 4, 6, 2,//------
594  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 1,//------
595  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2,//------
596  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 3,//------
597  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 4,//------
598  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 5,//------
599  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4,//------
600  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 1,//------
601  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 2,//------
602  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 1, 2,//------
603  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 2, 2,//------
604  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 3, 2,//------
605  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 2,//------
606  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 5, 2,//------
607  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 6, 2,//------
608  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 7, 2,//------
609  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 4, 2,//------
610  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 5, 2,//------
611  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 2,//------
612  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 2, 1,//------
613  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 2, 2,//------
614  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 2, 3,//------
615  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 2, 4,//------
616  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 2, 5,//------
617  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 2, 2, 4,//------
618  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 2, 2, 4, 1,//------
619  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 2, 2, 4, 2,//------
620  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 2, 2, 4, 1, 2,//------
621  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 2, 2, 2, 4, 2,//------
622  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 3, 2, 2, 4, 2,//------
623  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 4, 2, 2, 4, 2,//------
624  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 5, 2, 2, 4, 2,//------
625  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 2, 2, 4, 2,//------
626  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 7, 2, 2, 4, 2,//------
627  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 7, 2, 2, 4, 1, 2,//------
628  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 9, 2, 2, 4, 2,//------
629  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 10, 2, 2, 4, 2,//------
630  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 11, 2, 2, 4, 2,//------
631  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 12, 2, 2, 4, 2,//------
632  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 13, 2, 2, 4, 2,//------
633  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 2,//------
634  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 1, 2,//------
635  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 2, 2,//------
636  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 3, 2,//------
637  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 2,//------
638  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 5, 2,//------
639  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 6, 2,//------
640  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 7, 2,//------
641  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 9, 1,//------
642  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 1,//------
643  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 2,//------
644  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 2, 1,//------
645  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 2, 2,//------
646  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 2, 3,//------
647  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 2, 4,//------
648  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 2, 2, 3,//------
649  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 2, 2, 4,//------
650  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 2, 2, 4, 1,//------
651  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 2, 2, 4, 2,//------
652  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 2, 2, 4, 1, 2,//------
653  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 2, 2, 4, 2, 2,//------
654  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 2, 2, 2, 4, 1, 2,//------
655  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 3, 2, 2, 4, 1, 2,//------
656  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 4, 2, 2, 4, 1, 2,//------
657  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 6, 2, 2, 4, 2,//------
658  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 7, 2, 2, 4, 2,//------
659  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 7, 2, 2, 4, 1, 2,//------
660  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 8, 2, 2, 4, 1, 2,//------
661  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 10, 2, 2, 4, 2,//------
662  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 11, 2, 2, 4, 2,//------
663  2, 2, 2, 4, 2, 2, 4, 4, 6, 2, 2, 4, 4, 6, 6, 8, 2, 2, 4, 4, 6, 12, 2, 2, 4, 2 //-----
664 };
665 
667 {
668  0 , // nonexisting zero element
669 
670  1 , 1 , 2 , 2 , 3 , 3 , 4 , 4 , 3 , 4 , // 1 - 10
671 
672  5 , 5 , 6 , 6 , 6 , 6 , 6 , 7 , 8 , 8 , // 11 - 20
673 
674  9 , 9 , 9 , 9 , 9 , 9 , 9 , 10 , 10 , 10 , // 21 - 30
675 
676 11 , 11 , 11 , 11 , 11 , 12 , 13 , 13 , 14 , 14 , // 31 - 40
677 
678 14 , 14 , 14 , 14 , 14 , 15 , 15 , 15 , 16 , 16 , // 41 - 50
679 
680 // ----------------------------------------------------------
681 
682 16 , 16 , 16 , 17 , 18 , 18 , 19 , 19 , 19 , 19 , // 51 - 60
683 
684 19 , 19 , 19 , 20 , 19 , 19 , 19 , 19 , 19 , 20 , // 61 - 70
685 
686 21 , 21 , 21 , 21 , 21 , 21 , 21 , 21 , 22 , 22 , // 71 - 80
687 
688 23 , 23 , 23 , 23 , 24 , 24 , 25 , 25 , 26 , 26 , // 81 - 90
689 
690 27 , 27 , 27 , 26 , 26 , 27 , 27 , 26 , 26 , 26 // 91 - 100
691 
692 };