Ignore:
Timestamp:
Dec 16, 2009, 12:14:47 PM (16 years ago)
Author:
garnier
Message:

CVS update

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/documents/UserDoc/DocBookUsersGuides/ForApplicationDeveloper/xml/Fundamentals/biasing.xml

    r1211 r1222  
    22<!--                                                          -->
    33<!--  [History]                                               -->
     4<!--    New section on Reverse MC: L. Desorgher, Dec-2009     -->
    45<!--    Converted to DocBook: Katsuya Amako, Aug-2006         -->
    56<!--    Changed by: Dennis Wright, 25-Jun-2002                -->
    … …  
    10991100
    11001101</sect2>
     1102
     1103<!-- ******************* Section (Level#2) ****************** -->
     1104<sect2 id="sect.EvtBias.ReverseMC">
     1105<title>
     1106Adjoint/Reverse Monte Carlo
     1107</title>
     1108
     1109<para>
     1110Another powerful biasing technique available in Geant4 is the Reverse Monte Carlo (RMC) method,
     1111also known as the Adjoint Monte Carlo method. In this method  particles are generated  on the external
     1112boundary  of the sensitive  part of the  geometry and then are tracked backward in the geometry till they
     1113reach  the external source surface,
     1114or exceed an energy threshold. By this way the computing time is focused only on particle tracks
     1115that are contributing to the tallies.
     1116The RMC  method is much rapid than the Forward MC method when the sensitive part of the geometry is
     1117small compared to the rest of the geometry and to the external source,  that has to be  extensive and
     1118not beam like. At the moment the RMC method is implemented in Geant4 only for some electromagnetic
     1119processes (see <xref linkend="sect.EvtBias.ReverseMC.InG4.Process" />). An example illustrating the use
     1120of the Reverse MC method in Geant4 is distributed within the Geant4
     1121toolkit in <emphasis role="bold">examples/extended/biasing/ReverseMC01</emphasis>.
     1122</para>
     1123
     1124<!-- ******************* Section (Level#3) ****************** -->
     1125<sect3 id="sect.EvtBias.ReverseMC.InG4">
     1126<title>
     1127Treatment of the  Reverse MC method in Geant4
     1128</title>
     1129
     1130<para>
     1131Different G4Adjoint classes have been implemented into the Geant4
     1132toolkit in order to run an adjoint/reverse simulation in a Geant4 application.
     1133This implementation is illustrated in  <xref linkend="fig.EvtBias.ReverseMC.InG4_1" />.
     1134An adjoint run is divided in a serie of alternative adjoint and forward tracking  of adjoint and normal particles.
     1135One Geant4 event treats one of this tracking phase.
     1136</para>
     1137
     1138<figure id="fig.EvtBias.ReverseMC.InG4_1">
     1139<title>
     1140 Schematic view of an adjoint/reverse simulation in Geant4
     1141</title>
     1142<mediaobject>
     1143  <imageobject role="fo">
     1144    <imagedata fileref="./AllResources/Fundamentals/ReverseMC_tracking.png"
     1145               format="png"  align="center" />
     1146  </imageobject>
     1147  <imageobject role="html">
     1148    <imagedata fileref="./AllResources/Fundamentals/ReverseMC_tracking.png"
     1149               format="png"  align="center" />
     1150  </imageobject>
     1151</mediaobject>
     1152</figure>
     1153
     1154
     1155<!-- ******************* Section (Level#4) ****************** -->
     1156<sect4 id="sect.EvtBias.ReverseMC.InG4.AdjointTracking">
     1157<title>
     1158Adjoint tracking phase
     1159</title>
     1160
     1161<para>
     1162 Adjoint particles (adjoint_e-, adjoint_gamma,...) are generated  one by one on the so called adjoint source
     1163with  random position,  energy (1/E distribution) and direction. The adjoint source is the
     1164external surface of a user defined volume or of a user defined sphere. The adjoint
     1165source should contain one or several sensitive volumes and should be small compared to the entire geometry.
     1166The user can set the minimum and maximum energy of the adjoint source. After its
     1167generation the adjoint primary  particle is tracked backward in the geometry till a user  defined
     1168 external surface
     1169(spherical or boundary of a volume) or is killed before if it reaches a user defined  upper energy limit that
     1170represents the maximum energy of the external source. During the reverse tracking, reverse
     1171processes take place where  the adjoint particle being tracked can be either scattered
     1172or transformed in another type of adjoint particle. During the reverse tracking the
     1173G4AdjointSimulationManager replaces the user defined primary, run, stepping, ... actions, by its own actions.
     1174A reverse tracking phase corresponds to one Geant4 event.
     1175</para>
     1176
     1177</sect4>
     1178
     1179<!-- ******************* Section (Level#4) ****************** -->
     1180<sect4 id="sect.EvtBias.ReverseMC.InG4.ForwardTracking">
     1181<title>
     1182Forward tracking phase
     1183</title>
     1184
     1185<para>
     1186When an adjoint particle reaches the external surface its weight, type,  position,
     1187and direction are registered and a  normal primary particle, with a type equivalent to the last generated primary
     1188adjoint, is generated with the same energy, position but opposite direction  and is  tracked in the forward direction
     1189in the sensitive region as in a forward MC simulation.
     1190During this forward tracking phase the 
     1191event, stacking, stepping, tracking actions defined by the user for his forward simulation are used.
     1192By this clear separation between
     1193adjoint and forward tracking phases, the code of the user developed for a forward simulation
     1194should be only slightly  modified to adapt it for an adjoint simulation (see <xref linkend="sect.EvtBias.ReverseMC.UserCoding" />).
     1195Indeed  the computation of the signals is done by the same actions
     1196or classes that the one used in the forward simulation mode.   
     1197A forward tracking phase corresponds to one G4 event.
     1198</para>
     1199
     1200</sect4>
     1201
     1202<!-- ******************* Section (Level#4) ****************** -->
     1203<sect4 id="sect.EvtBias.ReverseMC.InG4.Process">
     1204<title>
     1205Reverse processes
     1206</title>
     1207
     1208<para>
     1209During the reverse tracking, reverse processes act on the adjoint particles.
     1210The reverse processes that are at the moment available in Geant4 are the:
     1211<itemizedlist spacing="compact">
     1212        <listitem><para>
     1213        Reverse discrete  ionization for e-, proton and ions
     1214        </para></listitem>
     1215        <listitem><para>
     1216        Continuous gain of energy by ionization and bremsstrahlung for e- and by ionization for protons and ions
     1217        </para></listitem>
     1218        <listitem><para>
     1219        Reverse discrete e- bremsstrahlung
     1220        </para></listitem>
     1221        <listitem><para>
     1222        Reverse photo-electric effect
     1223        </para></listitem>
     1224        <listitem><para>
     1225        Reverse Compton scattering
     1226        </para></listitem>
     1227        <listitem><para>
     1228        Approximated multiple scattering  (see comment in <xref linkend="sect.EvtBias.ReverseMC.Limitation.ReverseMS"/>)
     1229        </para></listitem>
     1230</itemizedlist>
     1231It is important to note that the electromagnetic reverse processes are cut dependent
     1232as their equivalent forward processes. The implementation of the reverse processes is based on
     1233the forward processes implemented in the G4 standard electromagnetic package.               
     1234</para>
     1235
     1236</sect4>
     1237
     1238<!-- ******************* Section (Level#4) ****************** -->
     1239<sect4 id="sect.EvtBias.ReverseMC.InG4.RemarkOnNbEvents">
     1240<title>
     1241Nb of adjoint particle types and nb of G4 events of an adjoint simulation
     1242</title>
     1243
     1244<para>
     1245The list of type of adjoint and forward particles that are generated on the adjoint source
     1246and considered in the simulation  is a function of the adjoint processes declared in the physics list.
     1247For example if only the e- and gamma electromagnetic processes are considered,
     1248only adjoint e- and adjoint gamma will be considered as primaries.
     1249In this case an adjoint event will be divided in four G4 event consisting in  the reverse tracking of an adjoint e-,
     1250the forward tracking of its equivalent forward e-, the reverse tracking of
     1251an adjoint gamma, and the forward tracking of its equivalent forward gamma.
     1252In this case a run of 100 adjoint events will consist into 400 Geant4 events.
     1253If the proton ionization is also considered adjoint and forward protons are also generated as primaries
     1254and 600 Geant4 events are processed for  100 adjoint events.
     1255</para> 
     1256
     1257</sect4>
     1258
     1259</sect3>
     1260
     1261<!-- ******************* Section (Level#3) ****************** -->
     1262<sect3 id="sect.EvtBias.ReverseMC.UserCoding">
     1263<title>
     1264How to update a G4 application  to use the reverse Monte Carlo mode
     1265</title>
     1266
     1267<para>
     1268Some modifications are needed  to an existing  Geant4 application in order to adapt
     1269it for the use of the reverse simulation mode
     1270(see also the G4 example <emphasis role="bold">examples/extended/biasing/ReverseMC01</emphasis>).
     1271It consists into the:
     1272
     1273<itemizedlist spacing="compact">
     1274        <listitem><para>
     1275        Creation of  the adjoint simulation manager in the main code
     1276        </para></listitem>
     1277        <listitem><para>
     1278        Optional declaration of user actions that will be used during the adjoint tracking phase
     1279        </para></listitem>
     1280        <listitem><para>
     1281        Use of a special physics lists that combine the adjoint and forward processes
     1282        </para></listitem>
     1283        <listitem><para>
     1284        Modification of the user analysis part of the code
     1285        </para></listitem>
     1286</itemizedlist>
     1287</para>
     1288
     1289
     1290<!-- ******************* Section (Level#4) ****************** -->
     1291<sect4 id="sect.EvtBias.ReverseMC.UserCoding.Main">
     1292<title>
     1293Creation of   G4AdjointSimManager in the main
     1294</title>
     1295<para>
     1296The class G4AdjointSimManager represents the manager  of an adjoint simulation.
     1297This static class should be created somewhere in the main code. The way to do that is illustrated below
     1298<informalexample>
     1299<programlisting>
     1300   int main(int argc,char** argv) {
     1301      ...
     1302      G4AdjointSimManager* theAdjointSimManager = G4AdjointSimManager::GetInstance();
     1303      ...
     1304   }
     1305</programlisting>
     1306</informalexample>
     1307By doing this  the G4 application  can be   run in the reverse  MC mode as well
     1308 as in the forward MC mode. 
     1309It is important to note that G4AdjointSimManager
     1310is not a new G4RunManager and that the creation of G4RunManager in the main  and the declaration of the geometry,
     1311 physics list,  and user actions to G4RunManager is still needed.
     1312The   definition of  the adjoint and external sources and the start of an adjoint simulation can be controlled
     1313by   G4UI  commands in the directory  <emphasis role="bold">/adjoint</emphasis>.
     1314</para>
     1315</sect4>
     1316
     1317
     1318<!-- ******************* Section (Level#4) ****************** -->
     1319<sect4 id="sect.EvtBias.ReverseMC.UserCoding.AdjointActions">
     1320<title>
     1321Optional declaration of adjoint user actions
     1322</title>
     1323<para>
     1324During an adjoint  simulation the user stepping, tracking, stacking and event  actions declared to G4RunManager
     1325 are used only during the G4 events dedicated to the forward tracking of normal particles in the sensitive region,
     1326 while during the events where  adjoint particles are tracked backward the   following happen concerning these actions:
     1327   
     1328<itemizedlist spacing="compact">
     1329<listitem><para>
     1330The user stepping action is replaced by G4AdjointSteppingAction that is reponsible to stop an adjoint track when
     1331it reaches the external source,  exceed the maximum energy of the external source, or cross the adjoint source surface.
     1332If needed the user can declare its own stepping action  that will be called by  G4AdjointSteppingAction after the
     1333check of stopping track conditions. This stepping action can be different that the stepping action used for the
     1334forward simulation. It is declared to G4AdjointSimManager by the following lines of code :
     1335       
     1336<informalexample><programlisting>
     1337   G4AdjointSimManager* theAdjointSimManager = G4AdjointSimManager::GetInstance();
     1338   theAdjointSimManager-&gt;SetAdjointSteppingAction(aUserDefinedSteppingAction);
     1339</programlisting></informalexample>
     1340                 
     1341</para></listitem>
     1342
     1343<listitem><para>
     1344No stacking, tracking and event actions are considered by default. If needed the user can declare to G4AdjointSimManager
     1345stacking, tracking and event actions that will be used only during the adjoint tracking phase.
     1346The following lines of code show how to declare these adjoint actions to G4AdjointSimManager:
     1347
     1348<informalexample><programlisting>
     1349   G4AdjointSimManager* theAdjointSimManager = G4AdjointSimManager::GetInstance();
     1350   theAdjointSimManager-&gt;SetAdjointEventAction(aUserDefinedEventAction);
     1351   theAdjointSimManager-&gt;SetAdjointStackingAction(aUserDefinedStackingAction);
     1352   theAdjointSimManager-&gt;SetAdjointTrackingAction(aUserDefinedTrackingAction);
     1353</programlisting></informalexample>
     1354
     1355</para></listitem>
     1356       
     1357</itemizedlist>
     1358By default no user run action is considered in an adjoint simulation but if needed such action can be declared to
     1359 G4AdjointSimManager as such:
     1360<informalexample><programlisting>
     1361   G4AdjointSimManager* theAdjointSimManager = G4AdjointSimManager::GetInstance();
     1362   theAdjointSimManager-&gt;SetAdjointRunAction(aUserDefinedRunAction);
     1363</programlisting></informalexample>
     1364
     1365</para>
     1366
     1367</sect4>
     1368
     1369<!-- ******************* Section (Level#4) ****************** -->
     1370<sect4 id="sect.EvtBias.ReverseMC.UserCoding.PhysicsList">
     1371<title>
     1372Physics list for reverse and forward electromagnetic processes
     1373</title>
     1374
     1375<para>
     1376To run an adjoint simulation a specific physics list should be used where existing G4 adjoint electromagnetic processes
     1377and their forward equivalent have to be declared.
     1378 An example of such physics list is provided by the class G4AdjointPhysicsLits in the  G4 example
     1379 <emphasis role="bold">extended/biasing/ReverseMC01</emphasis>.
     1380</para>
     1381
     1382</sect4>
     1383
     1384<!-- ******************* Section (Level#4) ****************** -->
     1385<sect4 id="sect.EvtBias.ReverseMC.UserCoding.Analysis">
     1386<title>
     1387Modification in the analysis part of the code
     1388</title>
     1389
     1390<para>
     1391The  user code should be modified to  normalize the signals computed during the forward tracking phase to the weight of the last adjoint particle
     1392that reaches the external surface.
     1393This weight represents the  statistical weight  that the last full adjoint tracks (from the adjoint source
     1394to the external source)  would have in a forward  simulation.   If multiplied by a signal  and registered in function of energy
     1395and/or direction the simulation results will give an answer matrix of this signal.
     1396 To normalize it to a given spectrum it has to be furthermore multiplied by a directional differential flux 
     1397 corresponding to this spectrum
     1398   
     1399The weight, direction, position , kinetic energy and type of the last adjoint particle that reaches the
     1400 external source, and that would represents the primary of a forward simulation,  can be   get from G4AdjointSimManager by using for example the following line of codes
     1401<informalexample><programlisting>
     1402   
     1403   G4AdjointSimManager* theAdjointSimManager = G4AdjointSimManager::GetInstance();
     1404   G4String particle_name = theAdjointSimManager-&gt;GetFwdParticleNameAtEndOfLastAdjointTrack();
     1405   G4int PDGEncoding= theAdjointSimManager-&gt;GetFwdParticlePDGEncodingAtEndOfLastAdjointTrack();
     1406   G4double weight = theAdjointSimManager-&gt;GetWeightAtEndOfLastAdjointTrack();
     1407   G4double Ekin  = theAdjointSimManager-&gt;GetEkinAtEndOfLastAdjointTrack();
     1408   G4double Ekin_per_nuc=theAdjointSimManager-&gt;GetEkinNucAtEndOfLastAdjointTrack(); // in case of ions
     1409   G4ThreeVector dir =  theAdjointSimManager-&gt;GetDirectionAtEndOfLastAdjointTrack();
     1410   G4ThreeVector pos = theAdjointSimManager-&gt;GetPositionAtEndOfLastAdjointTrack();
     1411   
     1412</programlisting></informalexample>
     1413
     1414</para>
     1415
     1416<para>
     1417In order to have a code working for both forward and adjoint simulation mode, the extra code needed in user actions or analysis
     1418manager  for the adjoint
     1419simulation mode can be separated to the code needed only for the normal forward  simulation by using the following public method
     1420of G4AdjointSimManager:
     1421<informalexample><programlisting>
     1422G4bool GetAdjointSimMode();
     1423</programlisting></informalexample>
     1424that returns true if   an adjoint simulation is running and false if not.
     1425
     1426</para>
     1427
     1428<para>
     1429The following   code example shows how  to normalize a detector signal and compute an answer matrix
     1430in the case of an adjoint simulation.
     1431 <example id="sect.EvtBias.ReverseMC.UserCoding.Analysis_1">
     1432<title>
     1433Normalization in the case of an adjoint simulation. The detector signal S computed during the forward tracking
     1434phase is normalized to a primary source of e-  with a differential directional flux given by the function F.
     1435An answer matrix of the signal is also computed. 
     1436</title>
     1437
     1438<programlisting>
     1439   
     1440   G4double S = ...; // signal   in the sensitive volume computed during a forward tracking phase
     1441
     1442   //Normalization of the signal for an adjoint simulation
     1443   G4AdjointSimManager* theAdjSimManager =    G4AdjointSimManager::GetInstance();
     1444   if (theAdjSimManager->GetAdjointSimMode()) {
     1445        G4double normalized_S=0.; //normalized to a given  e- primary spectrum
     1446        G4double S_for_answer_matrix=0.; //for e- answer matrix
     1447       
     1448        if (theAdjSimManager->GetFwdParticleNameAtEndOfLastAdjointTrack() == "e-"){
     1449            G4double ekin_prim = theAdjSimManager-&gt;GetEkinAtEndOfLastAdjointTrack();
     1450            G4ThreeVector dir_prim =  theAdjointSimManager-&gt;GetDirectionAtEndOfLastAdjointTrack();
     1451            G4double weight_prim = theAdjSimManager->GetWeightAtEndOfLastAdjointTrack();
     1452            S_for_answer_matrix = S*weight_prim;
     1453            normalized_S = S_for_answer_matrix*F(ekin_prim,dir); //F(ekin_prim,dir_prim) gives the differential directional flux of primary e-
     1454        } 
     1455       //follows the code where normalized_S  and S_for_answer_matrix are  registered or whatever
     1456        ....
     1457   }
     1458   
     1459   //analysis/normalization code for forward simulation
     1460   else {
     1461     ...
     1462   }
     1463   ...
     1464</programlisting>
     1465</example>
     1466
     1467</para>
     1468
     1469</sect4>
     1470
     1471</sect3>
     1472
     1473<!-- ******************* Section (Level#3) ****************** -->
     1474<sect3 id="sect.EvtBias.ReverseMC.Control">
     1475<title>
     1476Control of an adjoint simulation
     1477</title>
     1478
     1479<para>
     1480The G4UI  commands in the directory 
     1481<ulink url="./AllResources/Control/UIcommands/_adjoint_.html">/adjoint</ulink>.
     1482allow the user to :
     1483<itemizedlist spacing="compact">
     1484        <listitem><para>
     1485        Define the adjoint source where adjoint primaries are generated
     1486        </para></listitem>
     1487        <listitem><para>
     1488        Define the external source till which adjoint particles are  tracked
     1489        </para></listitem>
     1490        <listitem><para>
     1491        Start an adjoint simulation
     1492        </para></listitem>
     1493</itemizedlist>
     1494</para>
     1495</sect3>
     1496
     1497<!-- ******************* Section (Level#3) ****************** -->
     1498<sect3 id="sect.EvtBias.ReverseMC.Limitation">
     1499<title>
     1500Known issues in the Reverse MC mode 
     1501</title>
     1502
     1503<!-- ******************* Section (Level#4) ****************** -->
     1504<sect4 id="sect.EvtBias.ReverseMC.Limitation.Weight">
     1505<title>
     1506Occasional wrong high weight in the adjoint simulation
     1507</title>
     1508
     1509<para>
     1510In rare cases an adjoint track may get a wrong high weight when  reaching the external source.
     1511While this happens not often it may corrupt the simulation results significantly. This happens
     1512in some tracks where  both  reverse photo-electric and bremsstrahlung processes take place at low energy.
     1513We still need some investigations to remove this problem  at the level of physical adjoint/reverse processes.
     1514However  this problem can be solved  at the level of event actions or analysis in the user code by adding a test on the
     1515normalized signal during an adjoint simulation. An example of such test has been implemented in the Geant4 example
     1516<emphasis role="bold"> extended/biasing/ReverseMC01 </emphasis>.
     1517In this implementation  an event is rejected when the relative error of the computed  normalized energy deposited
     1518increases during one event by more than 50% while the computed precision  is already below 10%.
     1519</para>
     1520
     1521</sect4>
     1522
     1523<!-- ******************* Section (Level#4) ****************** -->
     1524<sect4 id="sect.EvtBias.ReverseMC.Limitation.ReverseBrem">
     1525<title>
     1526Reverse bremsstrahlung
     1527</title>
     1528
     1529<para>
     1530A difference between the differential cross sections used in the adjoint and forward bremsstrahlung
     1531 models is the source of a higher flux of >100 keV gamma in the reverse simulation compared to the forward simulation mode.
     1532In principle the adjoint processes/models should make use of the direct differential cross section to sample
     1533 the adjoint  secondaries and compute the adjoint cross section. However due to the way the effective differential cross section is
     1534 considered in the forward model G4eBremsstrahlungModel this was not possible to achieve for the reverse bremsstrahlung.
     1535Indeed the differential cross section used in G4AdjointeBremstrahlungModel  is obtained by the numerical derivation
     1536over the cut energy of the direct cross section provided  by G4eBremsstrahlungModel.
     1537This would be a correct procedure if the  distribution of secondary in   G4eBremsstrahlungModel 
     1538would match this differential cross section. Unfortunately it is not the case as independent  parameterization are used
     1539 in   G4eBremsstrahlungModel  for both the cross sections and the sampling of secondaries. (It means that in the forward case
     1540 if one would integrate the effective differential cross section considered  in the simulation we would not find back
     1541 the used cross section).
     1542 In the future we plan to correct this problem by using an extra weight correction factor after the occurrence of a reverse
     1543 bremsstrahlung. This weight factor should be the ratio between the differential CS used in the adjoint simulation and the
     1544one effectively used in the forward processes. As it is impossible to have a simple and direct  access to the forward differential CS
     1545 in G4eBremsstrahlungModel we  are investigating the feasibility to use  the differential CS considered in 
     1546 G4Penelope models.
     1547</para>
     1548
     1549</sect4>
     1550
     1551<!-- ******************* Section (Level#4) ****************** -->
     1552<sect4 id="sect.EvtBias.ReverseMC.Limitation.ReverseMS">
     1553<title>
     1554Reverse multiple scattering
     1555</title>
     1556 
     1557<para>
     1558For the reverse multiple scattering  the same model is used than  in the forward case.
     1559This approximation makes that the discrepancy between the adjoint and forward
     1560simulation cases can get to a level of ~ 10-15% relative differences in the test cases that we have considered.
     1561In the future we plan to improve   the adjoint multiple scattering models  by forcing the computation of
     1562 multiple scattering effect at the end of an adjoint step.
     1563</para>
     1564 
     1565</sect4>
     1566
     1567</sect3>
     1568</sect2>
    11011569</sect1>
Note: See TracChangeset for help on using the changeset viewer.