Home

xmds: e X tensible M ulti D imensional S imulator

1. FIGURE 10 4 Results for tla xmds The cross propagating components N and P may still be thought of as existing in the same space as the main vector component E and so they are declared as an extra vector in the lt field gt element However the equations governing the evolution of such cross vectors are not allowed to include any transverse derivatives i e they are not allowed to be PDEs Therefore through the main vector equations may be PDEs the cross propagating vector need not be transformed to Fourier space when the main vector is So although the cross vector components could be included as part of the main vector they are better defined as a separate vector for efficiency reasons In the SIIP integration algorithm the transverse evolution of the cross vector is calculated simultaneously with the forward evolution of the main vector but in all other algorithms the cross vector is calculated prior to calculating the main vector derivatives Thus the governing equations for the cross vector must be separated from those for the main vector This is done by including a lt cross_propagation gt element within the main lt integrate gt element and placing the cross vector equations within Also required within this element 162 MORE EXAMPLES is a list of the lt vectors gt that are to
2. FIGURE 3 1 Three dimensional plot in Matlab of a cosine modulated Gaussian pulse according to the advection equation Parameters used were v 1 0 0 1 zo 0 k 7 0 Notice that with the periodic boundary conditions that the pulse moves off one side of the figure and re enters from the opposite side 3 4 2 Scilab Using xsil2graphics we generate the Scilab script by the command xsil2graphics scilab advection xsil and then running the following commands in Scilab gt exec advection sci gt temp_di zeros 50 51 t 1 zeros 1 51 temp d2 zeros 50 51 x 1 zeros 1 50 amp 1 zeros 50 51 error amp 1 zeros 50 51 gt advectionl fscanfMat advection1 dat Error Info buffer is too small too many columns in your file gt temp_di advectioni 1 3 5 ADDING COMMAND LINE ARGUMENTS 63 temp d2 advection1 2 amp 1 advectioni 3 gt error_amp_1 advection 7 4 Se_il s tomp cll e 3 x 1 temas 8 5 il clear advectioni temp di temp d2 gt ploreca l p e ly amo il which should generate something similar to that in Figure 0 20 0 61 0 0 0 50 1 0 0 48 FIGURE 3 2 Three dimensional plot in Scilab of a cosine modulated Gaussian pulse according to the advection equation Parameters used were v 1
3. 342 Seilabl 4 4 4 3 5 Adding command line arguments 3 6 The xmds template option 4 Stochastic simulations and MPI 41 Without MP coros 4 HE A EDS EAE 4 1 1 Making the simulation and getting results 41 1 1 Matlab and Octavel 41 12 DScellab o e 23 cerrara ac eS 4 1 2 Making the simulation hard 42 With MEI s o A 9 3 29 A pe EUR X Ge RR 4 2 1 Example using LAM MPI IT Numerical Modelling Theory 5 Introduction 79 CONTENTS ix 6 Numerical Modelling Theory 81 AA II es aie ee Gace A 81 6 1 1 Boundary Conditions 44 4 44 43 sad 444 2 amp 4 82 Soe ER a ee ee oe eee 83 6 2 1 Ito vs Stratonovich Calculus do aed xXx Rex EOS EO 83 De PEG s Mie Mie 83 6 3 1 The Euler and Inverse Euler Methods 84 6 3 2 The Improved Euler Method 84 6 3 3 The Semi Implicit Method 85 6 3 4 The Fourth Order Runge Kutta Method 85 6 3 5 The Adaptive Fourth Order Runge Kutta Method 86 PME 87 6 4 1 Evaluating Transverse Derivatives 87 hse a TT 89 AA 89 90 plicit Picture 6 4 8 Interaction Picture Methods 99 6 4 9 The Semi Implicit method in the Interaction Picture 99 6 4 10 The Fourth Order Runge Kutta Method in the Interaction Pi
4. temperature ICD diffusion1 3 error temperature 1 diffusion1 4 gt gt w_1 2 ien cu Ci 3 3 x 1 temas 8 5 il clear diffusioni temp di temp d2 Plotting the data plot3d x 1 t 1 temperature 1 should give something along the lines of that in Figure 0 0 FIGURE 1 4 Three dimensional plot in Scilab of the diffusion of Gaussian pulse according to the diffusion equation Parameters used were amp 0 1 0 0 1 zo 0 Note that in both Figure 1 3 and Figure 1 4 we have an initial Gaussian pulse at t 0 which then spreads out and loses amplitude as t increases This is the expected evolution of the solution according to the diffusion equation The data points near the back of the graph close to t 1 are unlikely to be an accurate representation of the solution at this point and indeed are unlikely to be correct Nevertheless we have the expected behaviour and have now demonstrated sufficient xmds tags for you the user to be confident to go off and write your own simulations We wish you the very best of luck 32 STARTING FROM SCRATCH Extra and advanced features 2 1 Error checking The error checking feature of xmds is enabled by default and is controlled by using the lt error_check gt tag This is a boolean tag and expects either a yes or no entry In the context of xmds error checking means to run the simulation twice once through at the full s
5. 9 3 The field element field lt name gt main lt name gt dimensions t lt dimensions gt lt lattice gt 100 lt lattice gt domains 5 5 lt domains gt lt samples gt 1 lt samples gt Sq Possiigily more mes code gt lt A gt This element is compulsory as it is central to the problem definition it specifies the geometry of the field There are five assignments that xmds looks for in this element lt name gt lt dimensions gt lt lattice gt lt domains gt and lt samples gt The lt name gt assignment provides a name for the field It is not compulsory and will default to main It is envisaged that future versions xmds may allow sub fields within fields but at present only one lt field gt element is permitted The lt dimensions gt assignment lists the names of the transverse field dimensions From this point on the number of transverse dimensions is set to the number of names found in this assignment If this assignment is empty or absent then the field is assumed to be without transverse dimensions as in the kubo xmds example and xmds will not look for either 140 WORKED EXAMPLE NLSE XMDS a lt lattice gt assignment or a lt domains gt assignment The dimension names must be valid names as far as the C language is concerned see Section The lt lattice gt assignment is compulsory if there were one or more transverse dimension names found in the lt dimensions gt assignment
6. lt Global variables for the simulation gt lt globals gt lt CDATA const double v 1 0 const double x0 0 0 const double sigma 0 1 56 USING A TEMPLATE double k M_PI sigma 1 gt lt globals gt Notice that we ve not used the const keyword in front of the k declaration and as signment This is because we re deriving the value of the wave number from the standard deviation of the wave and so the value itself as far as C C is concerned isn t a constant If you re worried that this may be a problem further down the track don t because it isn t 3 3 3 The field to be integrated over The field as per normal has a name of main We have one dimension other than the propagation dimension and that is x which we want to put say 50 grid points i e lt lattice gt is set to 50 and as mentioned earlier we want the domain to go from 0 5 to 0 5 so the lt domains gt tag is set to 0 5 0 5 We want to sample the first point of this dimension so we set lt samples gt to 1 Within the lt vector gt assignments we only have one vector to define and since we have to call at least one vector main we ll use that It will have to have a complex type because we re going to be using Fourier space for part of the integration hence we set lt type gt to complex There is only one component we have to define and that is the field A that we ll be integrating over We therefore use A as the lt compon
7. 2 The moment group s dimensions must have the same name and be in the same order as those of the field 3 If a dimension is specified as being in z space k space in the moment group then it must be initialised in x space k space This can be done using the lt fourier_space gt tag 4 Each dimension of the initialisation moment group must have the same number of points as the corresponding dimension of the field and the start and end coordinates must be the same as those for the initialisation moment group In other words in strict mode the geometry of the initialisation moment group must be the same to within some small variation as that of the field Note that XSIL files generated by a lt breakpoint gt tag automatically satisfy conditions 1 and 2 if the dimensions of the two simulations are the same and in the same order In the loose geometry matching mode the last condition is relaxed to 4 The step size in each dimension in the initialisation moment group must be the same as the step size in the corresponding dimension of the field 2 6 INITIALISATION OF FIELD VECTORS FROM FILE 41 5 Some of the moment group grid points must overlap i e a vector with points at positions x 0 2 4 6 8 cannot be initialised from a moment group with points at positions z 1 3 5 7 9 Note that points that aren t initialised by the moment group are set to zero or can be initialised by the CDATA element if set The advantag
8. lt filename gt lt group gt Attributes optional format asciilbinary optional precision double single Subelement of Path to tag lt simulation gt lt output gt Description Container for elements describing what data should be output and how it should be output Accepts two optional attributes format and precision The format attribute defines the output format of the data The options available are ascii and binary with ascii being the default option The precision attribute defines the output data 216 LANGUAGE REFERENCE precision This can be either double or single with these options referring to double or single precision floating point numbers respectively This option is only meaningful when format is set to binary as the precision does not affect the output when ascii is chosen Example lt simulation gt lt output format binary precision single gt lt amas tags gt lt output gt lt simulation gt 11 23 1 filename output optional lt filename gt string lt filename gt Contains string Subelement of Path to tag simulation output filename Description Optional filename for output data Example lt simulation gt lt output gt lt filename gt nlse xsil lt filename gt lt output gt lt simulation gt 11 23 2 group required lt group gt xmds tags lt group gt Contains Subelement of Path to tag simulation
9. lt filename gt lt group gt 8 2 SYNTAX SUMMARY 131 lt group gt lt output gt yes lt sampling gt lt post_propagation gt lt sampling gt lt group gt yes lattice lt fourier_space gt vectors moments type and C code post propagation lt group gt no lt fourier_space gt lt moments gt and C Code Table 8 2 The structural elements Element name Used in Req May contain lt name gt lt simulation gt no string defaults to filename extn lt field gt no string defaults to main lt vector gt yes string lt arg gt yes string lt prop_dim gt lt simulation gt yes string lt cross_propagation gt lt error_check gt lt simulation gt no yes no defaults to yes lt stochastic gt lt simulation gt no yes no defaults to no lt use_mpi gt lt simulation gt no yes no defs to no lt MPI_Method gt lt simulation gt no Scheduling Uniform defs to Scheduling lt paths gt lt simulation gt integer read if stochastic lt seed gt lt simulation gt integer 2 reqd if stochastic lt noises gt lt simulation gt integer read if stochastic dimensions field no array of strings lattice field array of integers lt sampling gt same no as trans dims lt integrate gt yes integer lt domains gt lt field gt
10. Contains array of integers Subelement of Path to tag lt simulation gt lt output gt lt group gt lt sampling gt lt lattice gt Description A space separated list of integers each describing how many points to sample of each transverse dimension if greater than 1 see below One entry must exist for each transverse dimension If an entry is set to 0 then xmds integrates the moments over this dimension this will cause the output field to no longer be a function of this transverse dimension If an entry is set to 1 then xmds will sample the moments on a cross sectional slice of this dimension also causing the output field to lose this transverse dimension If this dimension is in normal space then xmds will extract the slice at the middle lattice point point number N 2 1 using integer division otherwise xmds will extract the slice at the zero momentum point k 0 Example lt simulation gt lt output gt lt group gt lt sampling gt lt lattice gt 50 lt lattice gt lt sampling gt lt group gt lt output gt lt simulation gt 220 LANGUAGE REFERENCE 11 23 2 1 5 moments sampling required moments string variableName string variableName lt moments gt Contains array of strings Subelement of Path to tag simulation output group sampling moments Description A list of strings of the names of the moments to sample These var
11. S yet more ends code to come gt lt simulation gt 1 1 3 The field element Now that we have set up the physical constants of our problem we need to describe the field that we are going to be integrating over in other words we need to specify a discretised version of z t y t and z t at the start of the problem To specify the field for the simple example we are studying here we only need to tell xmds the initial value of the field i e that 8 STARTING FROM SCRATCH when t 0 Note that in more complex situations to be studied later there are more tags to worry about The first tag is the lt field gt element This is a container for the other information that we are using to describe the field and is used very simply as follows lt f Te q gt lt U Wows zmde tags in here gt lt field gt We next give the name of the field this is supplied with the lt name gt tag note that this is a sub element of the lt field gt tag and so is different from the lt name gt tag of the lt simulation gt element If no name is given for the field it defaults to main however as mentioned before it is a good idea to specify a name for the field anyway We ll call it main to be consistent with other scripts you are likely to see and because this is the main field We next must tell xmds which of the field moments to sample directly after the field is initialised This is not obvious to those new to xmds howeve
12. We continue to use the Fourier Transform or Spectral method to determine the transverse derivatives The procedure is 1 Calculate if required 2 x 99 ah 3 ag a 4 p F 39 ki Z faal 5 aa N x aa p 6 x x ag a1 h 7T ay a 551a 8 p F L x9 ki Z ay 9 a N x ay p 10 x 22 az a2 h 11 a a bgiaq b3 2a 12 p F 2 0 k Z a 13 a N x a p 14 xe x a4 a3 h 6 4 NUMERICAL METHODS FOR PARTIAL DIFFERENTIAL EQUATIONS 97 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 ag a byjaq ba gay bisa p F C 2 k1 F ad aa N x aa P z a as a4 h a a 05184 b5 2a b5 39 bs 4Ag p F L 2 k1 F a ae N x ae P 6 z a ag as h a a bsia boa boa bo 4aq be sae p F L 2 k1 ag a N x ai p x x az ag h aj a bz 12 b7 22 b 39 b7 Ag b7 5a b sa p F L 2 k1 aj aj N x aj p z x ag az h a a bg 19 bg ga baza p F L 2 k1 Z a ay N x a p x x ag ag h ac a bg 1aa bo ai bo za boga p F L 2 k1 F a a N x ac P x x aio ag h aq a bios bio 6a 510 78 bio sas b10 9ac p F L
13. a cay 34 a a ay The f coefficients for the last function evaluation are derived from the original Cash Karp coefficients see table in the following way g C16 CiCa fi 1 B6a eq C1 ber ca 4 9 f2 be fs bes bsa c1 c3 dorlczCa C3C4 9 fa beich beach g fs bes c5 be1ca beac1 9 Je beaci beica 9 6 36 This method requires six extra copies of the field One ay for the initial field value three ax ay ar to act as a working fields while the original vector a and a collect the derivative contributions It also requires a similar size vector for the derivatives p 6 4 6 The Ninth Order Runge Kutta Method in the Explicit Pic ture We continue to use the Fourier Transform or Spectral method to determine the transverse derivatives The procedure is 1 Calculate if required 2 x So gah 3 a a 4 p F 5 k F as 5 ag N x aa p 6 x9 x ag a1 h 7 ay a b21aa 8 p F L 2 k1 F ay 9 ap N x ay p 10 22 22 a3 az h 11 a a b3 a b32a 94 NUMERICAL MODELLING THEORY 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 p F L 1 k F a a N X ac P x x a4 a3 h aq a 544a ba gay 543a p F C z ki F faal ay N x aa p x a as
14. lt output gt lt simulation gt 11 23 2 2 2 moments post_propagation required lt moments gt string variableName string variableName lt moments gt Contains array of strings Subelement of lt post_propagation gt Path to tag simulation output group post propagation moments Description The names of the moments with different names to the moments defined directly within the group element to be derived from the post processing Example lt simulation gt lt output gt lt group gt lt post_propagation gt lt moments gt pow_dens lt moments gt lt post_propagation gt lt group gt 222 LANGUAGE REFERENCE lt output gt lt simulation gt Part V Appendix 223 The XSIL input output syntax Following is an example of a multi dimensional multi component field written in XSIL format for more information see reference 2 The XSIL format is extensible and the format below has been specifically tailored for representing this type of field The parent element known as an lt XSIL gt data container contains one lt Param gt assignment and two lt Array gt assignments The one and only lt Param gt assignment n_independent is an integer and defines the number of independent dimensions of the field The first lt Array gt element variables defines the total number and names of all variables both independent and dependent The second
15. xmds e X tensible M ulti D imensional Simulator An open source numerical simulation package Version 1 6 5 P T Cochrane G Collecutt P D Drummond and J J Hope February 26 2008 i Abstract Writing codes for the simulation of complex phenomena is an art and science unto itself What with finding and using good algorithms actually writing the code debugging the code and testing the code not much time is left to actually investigate what it was you were initially out to look at This is where xmds comes in xmds allows you to write a high level description of the problem you are trying to solve usually a differential equation of some form it goes away and writes low level simulation code for you trying very hard to keep the code as efficient as possible compiles it and presents it ready to be run xmds is an acronym for eXtensible Multi Dimensional Simulator It is open source software released under the GNU General Public License and is written to assist in the calculation of a wide range of problems with application in physics mathematics and even finance and economics A high level description of the problem at hand is written in XML the extensible markup language and xmds transforms this into C language code This code can then be compiled by a C C compiler to produce a binary executable which solves the problem about as quickly and efficiently as might be achieved with code written by an expert xmds
16. 1 2 5 The output element a ow AA 26 126 DhedinabScHpRE 2s Gyo ss s s a th he oe s o w w N W S w w S SRY 27 92 39921 9 99 a 28 EXACT EU OS ee ies Brack ee vs NE A OR NOE A RON 29 t2 2 Dolab zx ewe hey he EO ACE kU w aC Uu d 30 2 Extra and advanced features 33 2 1 Error IA Foul Ve Y we C 33 E 34 2 3 BenchmarkiBl 24 5 4 4 a di SEES e dede s m deed ED 35 As AA EA A Rd 35 25 Binary III 37 vil viii CONTENTS 2 5 1 The format attribute 2 5 2 The precision attribute 4 240 93 axe ne 2 6 Initialisation of field vectors from file 2 6 1 Intialisation from an XSILflel 2 6 2 Input data layout for ASCII and binary formats 2 6 3 Importing complex data ve isis 4x 9o gu jue bree Ed W BORNE Bee O 28 Pr f rences 44 ds eb ela eo Due rs r s ee e ES 2 81 Turning preferences on and offf jc gb Eri de BORNE Bee TS 2 8 9 What ate the oplions a x eoo ox 9o Re ek es EDET OOOO vin el w aa eo eee E e eh 3 Using a template 3 1 The advection equation 2 0 2 0 dass de o Rn 3 2 The template Code e des e 99 die ve eS 3 3 Ripping it to bits ae a lt 0 4339 EX o Ry Xo ere 3 3 1 Global system parameters and functionality Hate ue o AAA IL ig Se PLU NET A a RU NEURAL ERST EUR ER BOE MS 3 4 Making the simulation and getting results 3 4 1 Matlab and Octavel
17. 1 is the fifth order and ar the fourth order solution The coefficients that appear are listed in table 6 1 To adjust the time step xmds calculates the relative difference between the two solutions A If partial differential equations are being solved it will determine the maximum of the relative errors of all grid points omitting points where the function is less than Ctnresh x Ma where Cthresn is a small threshold value see chapter 11 and M a is the peak value of the function across the grid If multiple functions are present each peak value is determined separately If the prescribed accuracy is Ag a time step of size 67 is accepted if Am lt Ag It is rejected if Am gt A4 and the step must be calculated again The size of the next time step Tn 1 is determined by A 15 0 9287 if Am lt Asa OTn41 6 27 Atol Am 0 9267 if Am gt Atol 6 4 NUMERICAL METHODS FOR PARTIAL DIFFERENTIAL EQUATIONS 87 1 Qj bi Ci cs 1 37 2825 378 27648 1 1 215135 d 3 3 3 9 250 18575 10 40 10 621 48384 4 3 3 _ 9 6 125 13525 5 10 10 5 594 55296 11 5 70 35 277 5 1 54 2 27 27 0 14336 6 7 1631 175 575 44275 253 512 8 55296 512 13824 110592 4096 1771 4 j i xe 3 4 5 Table 6 1 Cash Karp parameters for the embedded Runga Kutta method This algorithm will consume even more memory than the fourth order Runge Kutta but this easily pays off in situations
18. 4 5 in e There are of course a few caveats when performing such a mapping Finally the last boundary condition is the initial state of the components The values or form of all the components must be specified for some point in the dimension of propagation x 6 2 STOCHASTIC EQUATIONS 83 6 2 Stochastic Equations Stochastic differential equations also include noise sources x which may be real or com plex ma i ij j 5 59 09 F x alx g a x 6 x 6 7 These can make propagation a little more technically challenging but not impossible and enables a vast range of complex phenomena to be modelled As a result an individual propagation path becomes dependent on the underlying noise sources used and usually many such paths using independent noises are combined and averaged to obtain the desired result The noise sources are normally delta correlated in the propagation dimension though not always in the transverse dimensions and not always with each other so that generally Ex a9 z Cx x1 6 8 EE Ce 6 a 2 Do x u xi 6 9 where as usual the subscript denotes the transverse dimensions 6 2 1 Ito vs Stratonovich Calculus The integration of a stochastic equation as expressed above is not well defined in normal calculus In fact there are multiple ways in which such an integration can be defined as the limit of some discrete process and they are not equivalent The two main forms
19. 63 ene usb bee eee 2 71 Gah ont UAE ate pid a e eee 73 7 1 xmds main procedures 4 44 Bae ee Oe Ee dE a r a OE EO Reges 119 7 2 xmds class hierarchy s e ooo Ron RR RR mo ee e 120 8 1 xmds a functional diagram 4 4 aw dees eee ea o3 Rogo 3 126 9 1 Results for nlse xmds 000000008 0 138 10 1 Results for ndparamp xmds instar su REOS Row A 152 10 2 Results for kubo xmds 000000000 00004 156 10 3 Results for ibre xxmds 0000000000000 00 2 158 10 4 Results for tlaxmds 0000000000000 00004 161 xiii xiv LIST OF FIGURES List of Tables seso qe atl 87 s S S SSS es ee 125 8 2 The structural elements 131 gee es do ee e See eae M 132 od eM ee ee ee ee ee ae a 133 i So O Se S cy ay nee es 134 XV xvl LIST OF TABLES Part I Tutorial Starting from scratch It turns out that one of the most difficult things to do in xmds is to write a script from scratch which is why most users take either one of their own old scripts or borrow someone else s or one of the examples as a template However this doesn t mean to say that one is never going to have to write a script from scratch just that it s usually quite unlikely therefore we explain how to do this here If you re not interested and want to get coding more quickly then it s alright to skip to Chapter 3 and learn how to modify an existing script to your needs In this chap
20. CDATA chippy gphila gphila ie lt moment_group gt lt moment_group gt lt moments gt ippy ichippy lt moments gt lt integrate_dimension gt no no no lt integrate_dimension gt 10 5 HIGHDIM XMDS 165 lt CDATA ippy phila ichippy gphila 1 gt lt moment_group gt lt moment_group gt lt moments gt py ic lt moments gt lt integrate_dimension gt yes yes no lt integrate_dimension gt lt CDATA py phila ic gphila I lt moment_group gt lt moment_group gt lt moments gt ppy ichi lt moments gt lt integrate_dimension gt no yes yes lt integrate_dimension gt lt CDATA ppy phila ichi gphila Ia lt moment_group gt vectors main vc lt vectors gt lt CDATA const complex densi phiibxphiia const complex dens3 phi3b phi3a const double gdensi const double gdens3 gphila re gphila re gphila im gphila im gphi3a re gphi3a re gphi3a im gphi3a im dphita_dt L2p phita ixVir gamlossii gamloss132 2 gamloss12 2 dx dy dz phila i gameff gamlossi2 densi phiia ix gami3 gamloss132 dens3x phila ixphiib xphi3a i chippy ippy dphiib_dt L2n phiib i Vir gamlossil gamloss132 2 gamloss12 2 dx dy dz phiib i gameff gamlossi2 densi phiib i gam13 gamloss132 dens3 philb i phila phi3b dphi3a_dt L4p phi3a i V3r gamloss31 gamlossi32 2 gamloss32 2 dx dy dz phi3a i gam33 gamloss32 dens3 phi3a
21. and in what regimes the assumption that our window size is large enough that the implicit periodic boundary conditions are unimportant for a discussion of what this last part of the sentence means have a look at Section 2 8 Preferences As of xmds 1 3 1 there has been the ability to have preferences specific to a given user This gives the user more flexibility than before as they can control how simulations are built without having to recompile and reinstall xmds Also if xmds was installed by the root user a non root user of the system can modify how their simulations are built without having to reinstall the system binary and therefore alter how all other users programs are built 2 8 1 Turning preferences on and off By default preferences are on and are used if the preferences file can be found If the preferences file cannot be found or if use prefs gt is set to no then the default settings defined at configuration and installation of xmds will be used This is also the case for preference flags that are not specified in the preferences file the default settings will be used Therefore one doesn t need to specify all of the settings able to be used just the ones one wishes to change The lt use_prefs gt tag is set in the global configuration section of the xmds simulation script namely at the beginning with the lt name gt lt prop_dim gt etc tags If lt use_prefs gt is set to yes explicitly or not set at all then xmd
22. and xmds will look for the same number of positive integers here as there were transverse dimensions Each of these integers defines the number of points or steps to use in the corresponding transverse dimension and so xmds expects to find integers with a value of 2 or more in each case The lt domains gt assignment is similar except that it defines the domain range for each dimension This is done by entering a bracketed pair of real numbers for each dimension and xmds expects the range b a for any pair a b to be greater than 107199 Each domain range is divided according to Equation 9 2 in which n is the number of lattice points Note that the points stop one spacing short of the upper value this is due to the forced periodic boundary conditions ajal 1 TA AS niu La er 9 2 Finally the lt samples gt assignment tells xmds which moment groups Section to sample immediately after all the field s vectors have been initialised A sequence of 0 or 4 entries is expected as many as there are moment groups defined in the output xmds checks this once all moment groups have been processed The correspondence is sequential A 1 means sample the moment group in question a 0 means do not Next xmds will look for the field s vectors 9 4 The vector element lt vector gt lt name gt main lt name gt lt type gt complex lt type gt lt components gt phi lt components gt lt fourier_space gt no l
23. b 59 b7 ea 1 ar hL 2 k a e aj aj E IN x F7 a6 106 NUMERICAL MODELLING THEORY 36 ay et ahe ose aj 37 x x ag a7 h 38 ay a bg 9 bg ga bs7a Sap e las h 2 k1 E 40 ay h JF N 22 FA aj dp eU a8 h x k1 TT 42 1 x ag ag h 43 a a by 1aa boa 097a bo gay ase e Uao hL 1 k EN 45 ace FIN OO m 46 a eti as h ak ac 47 x x aio ag h 48 aq a biota bio sa b10 7a bio sa b10 9a 49 Ad e Uaro hL 2 k Fa 50 aa h F INT FF faa 63 51 Ad eU amp o h x k1 T 52 x x ay a19 h 53 ae a b11 19 bita b11 729 b11 sas b11 99 01131084 54 a e C IDORE eL a 55 ae h F N x F7 ae E as eO mh a e 57 x x aiz a11 h 98 ay a bio 192 bio sa 01578 bio say bio 92 0151984 bio 119 59 af e Uar2 hL 20 k ay 60 ag h F N 2 las E 6 4 NUMERICAL METHODS FOR PARTIAL DIFFERENTIAL EQUATIONS 107 61 ay eel a 62 x x a3 ai2 h 63 a a b131a bisa 013 78 bis sas bis 99 bis 104 013 186 bis 1994 64 a e U ms3 h x0 k1 a 65 a h F INF lla 0 66 a eO ms a9 ay 67 x x a4 a13 h 68 a a b14 1aa b14 69 b14 7a bia say 014 986 D14 108q 4 014 1186 b14 1297 bisa 69 ay e monte a
24. bool lt use_openmp gt Contains boolean Subelement of Path to tag lt simulation gt lt use_openmp gt Description Defaults to no Tells xmds whether to make use of compiler di rectives to instruct the compiler how to parallelise parts of the simulation other than the FFT s The simulation will compile correctly if OpenMP by the compiler and you will need to ensure that the correct flags are passed to the compiler to ensure that it does enable OpenMP support e g for Intel s C compiler the flag openmp must be passed to the compiler to enable OpenMP OpenMP cannot be used simultaneously with MPI for deterministic simulations Also the number of OpenMP threads used in a simulation defaults to the number of physical processors available Note that best performance will be achieved if FFTW is compiled to use OpenMP threads instead of the default pthreads if your simulations make use of OpenMP threads See the OpenMP website http www openmp org for more information about OpenMP As of mid 2006 Intel s icc compiler supports OpenMP however the current release of GCC 4 1 1 does not Support for OpenMP is planned for the 4 2 release of GCC Example 186 LANGUAGE REFERENCE lt simulation gt lt use_openmp gt yes lt use_openmp gt lt simulation gt 11 18 fftw version optional fftw version int fftw version Contains integer Subelement of Path to tag simulation fftw version Des
25. gt gt kubo tutorial gt gt plot t 1 mean realz 1 gt gt xlabel t gt gt ylabel z generates Figure 0 9 4 0 5 4 04 J 0 3 4 0 2 l l 1 1 l 1 l l 1 0 FIGURE 4 1 Matlab generated plot of the evolution of the oscillator frequency z over time perturbed by noise 4 1 1 2 Scilab Using xsil2graphics we generate the Scilab script by the command xsil2graphics scilab kubo tutorial xsil Running the following sequence of commands in Scilab exec kubo tutorial sci gt temp_d1 zeros 1 101 t 1 zeros 1 101 mean realz 1 zeros 1 101 72 STOCHASTIC SIMULATIONS AND MPI gt mean_imagz_1 zeros 1 101 gt sd_realz_1 zeros 1 101 gt sd_imagz_1 zeros 1 101 error realz 1 zeros 1 101 error imagz 1 zeros 1 101 gt kubo o e A II E ASA A o d si ie dates Error Info buffer is too small too many columns in your file gt temp_di kubo_tutorial1 1 mean realz 1 kubo_tutoriali 2 VY Il mean imagz 1 kubo_tutorial1 3 sd realz 1 kubo_tutorial1 4 sd imagz 1 kubo tutoriali 5 error realz 1 kubo_tutorial1 6 error imagz 1 kubo_tutoriali 7 nt lC tempodi lt clear kubo_tutoriali temp_di generates Figure Notice that the results we get here are similar to the other Kubo oscillator example mentioned both later in this document and
26. i lt vector gt lt vector gt lt name gt main lt name gt 164 MORE EXAMPLES lt type gt complex lt type gt lt components gt phila phiib phi3a phi3b gphila gphi3a lt components gt lt fourier_space gt no no no lt fourier_space gt vectors vcl lt vectors gt lt CDATA const double realfn mu 0 5x Mx vcore Uoh11 hbar phila realfn gt 0 complex sqrt realfn 0 complex 0 0 phiib realfn gt 0 complex sqrt realfn 0 complex 0 0 phi3a complex 0 0 phi3b complex 0 0 gphiia realfn gt 0 complex sqrt realfn 0 complex 0 0 gphi3a complex 0 0 is lt vector gt lt field gt The sequence of integrations to perform gt lt sequence gt lt integrate gt lt algorithm gt ARK89IP lt algorithm gt lt interval gt 1e 7 lt interval gt lt tolerance gt 1 0e 7 lt tolerance gt lt lattice gt 1000 lt lattice gt lt samples gt 10 10 1 lt samples gt lt k_operators gt lt constant gt yes lt constant gt lt operator_names gt L2p L2n L4p L4n lt operator_names gt lt CDATA L2p complex 0 hbar M 2 chix kx kx ky ky kz kz L2n complex 0 hbar M 2 chix kx kx ky ky kz kz L4p complex 0 hbar M 4 chix kx kx ky ky kz kz L4n complex 0 hbar M 4 chix kx kx ky ky kz kz NES X k operators moment group lt moments gt chippy lt moments gt lt integrate_dimension gt yes yes yes lt integrate_dimension gt lt
27. i 0 5 phila xphila i gami3 gamlossi132 densi phi3a dphi3b dt L4n phi3b i V3r gamloss31 166 MORE EXAMPLES gamloss132 2 gamloss32 2 dx dy dz phi3b i gam33 gamloss32 dens3 phi3b i 0 5 phiib phiib ixgam13 gamloss132 densi phi3b dgphila_dt L2p gphiia ixgVir gamlossi1 gphiia t i gameff gamlossi2 gdensi ichippy i gami3 gamlossi32 gdens3 gphilia i conj gphiia gphi3a dgphi3a dt L4p gphi3a i gV3r gamloss31 gphi3a i gam33 gamloss32 gdens3 gphi3a i 0 5 gphilaxgphila t i gami3 gamlossi32 gdensi gphi3a gt lt integrate gt lt sequence gt lt I The output to generate gt gt gt lt output format binary precision double gt lt group gt lt sampling gt Stoner space no no noX fourier space lattice 16 1 1 lt lattice gt lt moments gt atoms molecules gatoms gmolecules lt moments gt lt CDATA atoms philbx phila molecules phi3b phi3a gatoms conj gphila gphila gmolecules conj gphi3a gphi3a 11 gt lt sampling gt lt group gt lt group gt lt sampling gt lt fourier_space gt no no no lt fourier_space gt lt lattice gt 0 16 0 lt lattice gt lt moments gt rn_1 rn_2 grn_1 grn 2 excitedn lt moments gt lt CDATA iu phiib xphila rn_2 phi3b phi3a grn_1 conj gphila gphiia grn_2 conj gphi3a gphi3a excitedn g g 4 phiib xphilb phila philat F F phi3b phi3a F g 2 phiib phiib phi3a phiia phiia phi3b 11 gt lt
28. lt Array gt element data defines the lattice for the data and then the data itself Following standard C index convention the last dimension index is the one that changes most rapidly in the data Hence this data has 3 lattice points in the first independent dimension t 5 lattice points in the second independent dimension x and there are 4 variables to be read 2 independent and 2 dependent which of course must be the same as the number of variable names declared in the variables array lt xml version 1 0 gt lt XSIL Name moment_group_1 gt lt Param Name n_independent gt 2 lt Param gt lt Array Name variables Type Text gt lt Dim gt 4 lt Dim gt lt Stream gt lt Metalink Format Text Delimiter n gt t x real al imag al lt Stream gt lt Array gt Array Name data Type double gt lt Dim gt 3 lt Dim gt lt Dim gt 5 lt Dim gt 225 226 THE XSIL INPUT OUTPUT SYNTAX lt Dim gt 4 lt Dim gt lt Stream gt lt Metalink Format Text Delimiter Nn gt 1 000000e 00 1 000000e 00 0 000000e 00 0 000000e 00 1 000000e 00 0 500000e 00 0 000000e 00 0 000000e 00 1 000000e 00 0 000000e 00 0 000000e 00 0 000000e 00 1 000000e 00 0 500000e 00 0 000000e 00 0 000000e 00 1 000000e 00 1 000000e 00 0 000000e 00 0 000000e 00 0 000000e 00 1 000000e 00 0 000000e 00 0 000000e 00 0 000000e 00 0 500000e 00 1 000000e 00 1 000000e 00 0 000000e 00 0 000000e 00 1 000000e 00 1 00000
29. or failing the existence of such a directory in the directory local to the simulation binary executable for a file labelled hostname wisdom where hostname is the name of the computer currently running the simulation If such a file exists then the xmds simulation will load this wisdom as part of the fftw plan creation step thus drastically reducing startup time of the simulation Any accumulated wisdom will be then added back to this file at the fftw plan deletion stage If no wisdom file is found the simulation will create one in an appropriate location and then save any accumulated wisdom to this file Example lt simulation gt lt use_wisdom gt no lt use_wisdom gt lt simulation gt 11 14 use double optional lt use_double gt bool lt use_double gt xmds 1 2 OBSOLETE see lt output gt Contains boolean 184 LANGUAGE REFERENCE Subelement of Path to tag simulation use double Description Defaults to yes Decides the precision of the output data Only useful when use binary is set to yes If set to no then single precision output is used This option is useful in reducing the size of the output files This tag is now obsolete Please specify the precision attribute in the output tag For example output precision double gt Example lt simulation gt lt use_double gt no lt use_double gt lt simulation gt 11 15 use prefs optional lt use_prefs gt bool lt us
30. p FL 2 k1 Z al b a AN x a p c a a ja p PA ki Z al a hN x a p ra a aj a 19 z h Here an extra copy of the field is needed to store the initial a into ar and a similar size vector is needed for the derivatives p 6 4 4 The Fourth Order Runge Kutta Method in the Explicit Pic ture Again we use the Fourier Transform or Spectral method to determine the transverse deriva tives The procedure after Section 6 3 4 then is 1 2 3 10 11 Calculate if required ax a a a po EG as ax h N x ax p a a ciax z ih ay a jay p F Ea ki ax ax hN x ax p a a sax 6 4 NUMERICAL METHODS FOR PARTIAL DIFFERENTIAL EQUATIONS 91 12 ax a lay 13 p F L 39 k1 F ax 14 ax hN x ag p 15 a a lag 16 x 2 ih 17 a ar a 18 p F L x k1 F ax 19 ax hN x ax p 20 a a iak This method requires two extra copies of the field One ay for the initial field value and the other ax to act as a working field while the original vector a collects the derivative contributions This saves the final step of copying ax back into a that would had to have been performed if ax had collected the derivative contributions It also requires a similar size vector for the derivatives p On most systems where the calculation is not memory limited this method tends to achieve a des
31. sampling gt lt group gt lt group gt lt sampling gt 10 6 HIGHDIM_VECTOR_VERSION XMDS 167 S four Tens pace no no no lt fourier_space gt lt lattice gt 4 8 16 lt lattice gt lt moments gt atomsr moleculesr atomsi moleculesi lt moments gt lt CDATA atomsr phila moleculesr phi3a atomsi ix gphila moleculesi i gphi3a 1 gt lt sampling gt lt group gt lt output gt lt simulation gt This simulation is included to highlight the usage of moment groups in evolution Here we wish to use various variables integrated over one two or all the transverse dimensions This is done in integrate or filter elements by the inclusion of a moment group Using the first as an example lt moment_group gt lt moments gt chippy lt moments gt lt integrate_dimension gt yes yes yes lt integrate_dimension gt lt CDATA chippy gphilaxgphila 11 gt lt moment_group gt The syntax is similar to output and filter syntax but note that the equality must be a and not simply This element provides the integral of the modulus squared of the field gphila which can then be used normally in the integration code by the designated name chippy Other variables are actually fields in which any number of transverse dimensions may be integrated and the others are left When lt functions gt and lt moment_group gt elements are used the position of the lt vectors gt tag is crucial It specifie
32. tag to yes We do this because if xmds has to assume that the k operators aren t constant then it has to use much slower code to evolve the solution and we are fortu nate that for our simulation here the k operators are constant since they can be calculated via the also constant x variable xmds needs to know the name of the operator we are going to use for our k operator and this we set with the lt operator_names gt tag to be L In general the lt operator_names gt tag expects a space separated list of the operator names you wish to define We then give the C code necessary to define the operator in a CDATA block and for our simulation this is lt CDATA L kappa kx kx 11 Note that we have one variable here that we haven t defined before kx This is a variable automatically defined by xmds when we define k operators If we had another variable called y and defined a k operator for it then it would be called ky The last thing we need to do within this section of the script is tell xmds how to evolve the solution in other words the differential equation As in Section we use a CDATA block with a modified C syntax that xmds understands to write down the differential equation lt CDATA dT_dt L T 11 gt This completes the work necessary to integrate the solution forward and completes the lt integrate gt and lt sequence gt elements which are The sequence of integrations to perform gt lt
33. 0 0 1 zo 0 k m c Notice that with the periodic boundary conditions that the pulse moves off one side of the figure and re enters from the opposite side 3 5 Adding command line arguments Now that we re happy with how the simulation is performing we can define some command line arguments so that we can investigate the system more easily We might as well be able to change all of the global variables except for k since it depends on one of the other variables so add the following code snippet just before the lt globals gt block 64 USING A TEMPLATE AI Command Whine angumenitis gt lt argv gt lt arg gt lt name gt v lt name gt lt type gt double lt type gt lt default_value gt 1 0 lt default_value gt lt arg gt lt arg gt lt name gt x0 lt name gt lt type gt double lt type gt lt default_value gt 0 0 lt default_value gt lt arg gt lt arg gt lt name gt sigma lt name gt lt type gt double lt type gt lt default_value gt 0 1 lt default_value gt lt arg gt lt argv gt But remember to comment out these variables in the lt globals gt block otherwise the compiler will throw an error Now you can have a play Try different variable options and see how the equations and xmds perform 3 6 The xmds template option A new feature of xmds as of version 1 3 3 is the output of a template code When xmds is called with the template or equivalently the t option t
34. 3 in which the argument of the complex vector z is blown about by a real Gaussian noise term t This is a simple stochastic ODE i t z 10 3 Such Gaussian noise terms in analytic form are correlated in time and space through Dirac delta functions as shown in Equation 10 4 amp x 6 x HT 08 x x 10 4 However when solving stochastic DEs numerically algorithms work with discrete time in tervals and lattice spacings Therefore these Dirac delta correlations must be transformed to Kronecker delta correlations using the integration time step and the spatial volume of the 154 MORE EXAMPLES lattice as shown in Equation 10 5 i IIS oua amp x 6 x TD Ant 10 5 The good news is that xmds calculates this for the user all that has to be done is to specify as a simple lt noises gt assignment within the lt simulation gt element the maximum number of noise terms required in any one segment and then reference them as n_1 n_2 etc as can be seen in this example These noises are available within the initialisation code for each field lt vector gt within the main integration equations code not in the lt k_operators gt code and in the code for any lt filter gt segments Within the initialisation code and lt filter gt code the variances of the noises are determined by the lattice cell volume product for the particular lt fourier_space gt specification as shown in Equ
35. 4 t h This is known at the explicit or forward time Euler method It is very simple but not very accurate or stable Stability can be improved to an extent by thinking in reverse Rather than using the slope at the current point in time to calculate the next value of a use the slope and the next point in time But this of course is not known until the next value of a is known so the process becomes an iterative or implicit one Initially the current value of a is used as next value and then the step Angi An h f ba ana 6 14 is iterated until convergence is obtained This method is known as the implicit or backward time or inverse Euler method For certain classes of differential equation implicit methods are more stable than their explicit equivalents but overall the implicit Euler method is no more accurate than the explicit Euler method The explicit and implicit Euler methods do not account for the second and higher order derivatives of a Thus the error per step is of order h which means the error over a given integration interval is of order h Consequently these methods are known as first order methods Both of these methods can be used to integrate stochastic equations by treating the noise terms as part of a larger function f which depends on a pseudo random number which changes each time step This technique converges to the Ito integral with a convergence of order h1 In general the order of the stochastic c
36. a binary file for input to xmds on PowerPC will not load correctly on an x86 platform For the binary and ascii formats xmds expects the input data to be essentially interlaced so that it is loaded into memory correctly The best way to explain this is by way of example If the input data is to be three vectors say x y and z then xmds expects the data to be formed thus x 0 y 0 z 0 x 1 y 1 z 1 and so on with a new line or space between each entry The only difference between the ascii and binary input formats is that the newline or space character is unnecessary and xmds just expects a sequence of double values ordered as above 11 21 FIELD 195 For complex types xmds expects the imaginary part of each variable to follow the real part If the input data for the previous example were complex numbers then xmds would expect the data to be in the order real x 0 imag x 0 real y 0 imag y 0 real z 0 imag z 0 Example lt simulation gt lt field gt lt vector gt lt filename format binary gt blah lt filename gt lt vector gt lt field gt lt simulation gt 11 21 6 3 type vector optional lt type gt string of complex or double lt type gt Contains string interpreted as either complex or double type Subelement of Path to tag simulation field lt vector gt type Description The data type of the vector Defaults to complex if not specified It is a good idea
37. and hence the corresponding xmds code is similarly nested The syntax of adding command line arguments to simulations is as follows lt argv gt lt arg gt lt name gt lt name gt lt gt the argument name gt lt type gt lt type gt lt data type of arg gt lt default_value gt lt default_value gt lt the default value gt lt arg gt more arg definitions here if necessary gt lt argv gt We ll now go through an example to show you how to use this feature and some of the subtleties of using command line arguments with xmds derived simulations Let s revisit the diffusion simulation discussed in Chapter 1 Section The main use of command line arguments is to be able to replace variables given in the lt globals gt element Therefore let s change the diffusion coefficient kappa in the code to be an argument to the simulation We do this by adding the following code before the lt globals gt element and by commenting out the kappa declaration and assignment in the existing code The xmds code then becomes simulation I global parameters and functionality tags in here gt lt I Command line arguments lt argv gt lt arg gt 2 7 COMMAND LINE ARGUMENTS 45 lt name gt kappa lt name gt lt type gt double lt type gt lt default_value gt 0 1 lt default_value gt lt arg gt lt argv gt lt l Global var ables for the simulation gt gt gt lt globals
38. array of bracketed pairs of floats lt components gt lt vector gt yes array of strings lt type gt lt vector gt no lt complex gt or lt double gt lt sampling gt lt arg gt yes string lt vectors gt lt vector gt no array of strings lt k_operators gt lt sampling gt lt integrate gt yes array of strings lt filter gt lt cross_propagation gt lt moment_group gt lt integrate gt no array of strings lt filter gt lt functions gt lt integrate gt no array of strings lt filter gt lt constant gt lt k_operators gt no yes no defaults to no 132 FUNCTIONALITY lt operator_names gt lt k_operators gt yes array of strings lt filename gt lt vector gt no string lt output gt no string defaults to lt name gt xsil lt fourier_space gt lt vector gt lt filter gt no array of yes no lt sampling gt yes array of yes no lt post_propagation gt yes array of yes no same no as non collapsed trans dims 1 lt samples gt lt field gt yes integer as many as are lt group gt lt integrate gt lt cycles gt lt sequence gt yes integer not in main sequence lt algorithm gt lt integrate gt no string defaults to RK4EX or SIEX for stochastic problems lt interval gt lt integrate gt yes string lt no noise gt lt integrate gt no yes no defaults to no lt iterations
39. be cross propagated and a lt prop_dim gt assignment specifying the dimension of cross propagation The lt vectors gt that were made accessible for the main equations will also be accessible here 10 5 highdim xmds Sram vensi on il OY eS lt simulation gt lt name gt highdim lt name gt gt gt Gioball system spananeners ands fine tionality lt prop_dim gt t lt prop_dim gt lt error_check gt yes lt error_check gt lt use_mpi gt yes lt use_mpi gt lt use_wisdom gt yes lt use_wisdom gt lt benchmark gt yes lt benchmark gt lt Global variables for the simulation gt lt globals gt lt CDATA const double noise 0 0 const double hbar 1 05500000000e 34 const double M 1 409539200000000e 25 const double omegax 0 58976353090742 const double omegay 0 58976353090742 const double omegaz 0 58976353090742 30 const double U11 2 974797272874263e 51 const double U13 1 417820412490823e 50 const double U33 2 974797272874263e 51 const double inum 1 0e6 const double Uohii Uii hbar const double Uohi3 U13 hbar const double Uoh33 U33 hbar const double mu pow 15 inum U11 omegax omegay omegaz M PI 4 0 4 pow M 0 6 2 const double delta 1 0e9 const double F 2 0e 2 const double g sqrt Uohi1 2 0 delta const double lossii 1 0e 2 const double lossi2 1 6e 22 const double loss31 1 0e 2 const double loss32 1 6e 22 const double loss1i32 8 0e 17 const double
40. boolean lt fourier_space gt Contains array of booleans Subelement of Path to tag lt simulation gt lt integrate gt lt filter gt lt fourier_space gt Description Tells xmds in which space the filter is to be applied This is a list of yes no options for each vector Example 212 LANGUAGE REFERENCE lt simulation gt lt integrate gt lt filter gt lt fourier_space gt yes no lt fourier_space gt lt filter gt lt integrate gt lt simulation gt 11 22 3 5 cross propagation optional lt cross_propagation gt xmds tags lt cross_propagation gt Contains prop CDATA Subelement of Path to tag lt simulation gt lt integrate gt lt cross_propagation gt Description Container of the tags describing the cross propagation vectors if any Example lt simulation gt lt integrate gt lt cross_propagation gt Sl sonde tags gt lt cross_propagation gt lt integrate gt lt simulation gt 11 22 3 5 1 prop dim cross propagation required lt prop_dim gt string lt prop_dim gt Contains string Subelement of Path to tag simulation integrate cross propagation prop dim Description The propagation dimension of the cross propagating vector Example lt simulation gt lt integrate gt lt cross_propagation gt lt prop_dim gt z lt prop_dim gt 11 22 SEQUENCE 213 lt cross_propagation gt lt integrate gt lt simulation
41. changes and one addition the lt fourier space gt tag The vector lt name gt is again main the lt type gt this time is complex though The reason this might be confusing is because the temperature is a real quantity and therefore those who have been reading carefully may question why float wasn t chosen as the type instead A type of complex is chosen because we are using a Fourier transform technique whereby the solution is transformed into Fourier space to calculate the second partial derivative in x and then transformed back again and to be able to transform into Fourier space we need our variables to be of complex type The lt components gt tag is set to T since this is the name of the variable about to be defined in the CDATA block to come and we tell xmds that this component is not defined in Fourier space by setting the lt fourier_space gt tag to no We next define the CDATA block This is our C language representation of Equa tion 1 6 which we are using as our initial condition The CDATA block is lt CDATA T rcomplex exp x x0 x x0 2 0 sigma sigma sigma sqrt 2 0 M_PI 70 0 gt This may look quite complicated and possibly because we ve tried to split the code over several lines in an attempt to break up the various parts and because this is intended for a fixed width page The code does nevertheless introduce some important concepts These are the rcomplex function the ability to
42. double array xOut_1 1x201 1608 double array yOut 1 1x201 1608 double array zDut_1 1x201 1608 double array Grand total is 1407 elements using 11256 bytes We can see that we ve loaded our output data from the simulation into the variables x0ut_1 yOut 1 and z0ut_1 in Matlab or Octave The reason why the _1 is appended to the variable names is so that if one defines a variable in two different moment groups but of the same name then data isn t lost The number refers to the label of the moment group in the simulation At present you don t need to worry about these details but just realise that they are there for when you write more complex scripts in the future The error variables seen in the whos listing are the differences between the full step integration and the half step integration The half step integration is used for error checking purposes so that you can check if your simulation is likely to be giving you reliable answers Now plot the data by going gt gt plop ne 1 wOwe_il wwe s gt gt xlabel x gt gt ylabel y gt gt z apel 22 and you should see a figure similar to Figure FIGURE 1 1 Three dimensional plot in Matlab of the trajectories of a Lorenz attractor Pa rameters used were c 10 b 8 3 r 28 with initial conditions of zo 1 0 yo 1 0 and 20 20 0 1 1 A BASIC SIMULATION 19 1 1 7 2 Scilab A very similar process is necessary for viewing the results in Scilab Star
43. enable fftw3 to check for version 3 of fftw and enable use of it for building simu lations e with fftw path the path to the fftw installation Note when fftw is installed in usr local then with fftw path must be specified so that the configure script can find the libraries and headers e with fftw3 path the path to the fftw3 installation if different to the path for fftw2 e with mpi libs extra libraries needed when checking for MPI e with mpi path set the path to the prefix of your MPI distribution e with mpi compiler set the mpi C compiler 8 1 2 Usage For usage at the command prompt type xmds and the output should be This is xmds version 1 6 5 using C compiler gcc and C compiler mpicc for parallel work Usage xmds v c infile infile required the input file v optional verbose mode c optional turns off automatic compilation of simulation 8 1 INSTALLING AND RUNNING xmds 129 or something similar The verbose mode can be quite useful as it will repeat back to the user exactly how it is interpreting the input file Switching off automatic compilation can be handy when one is testing a script and doesn t want to wait for the script to compile especially if the resultant binary isn t wanted anyway This option can also be used by those people who want to use xmds to produce a skeleton C code with which they want to change directly themselves to perform something m
44. equations What it does is it looks for operator component combinations registers their presence and the replaces their text with something else In the explicit picture they are replaced with a reference to an internally calculated vector p F L x9 k Z a but in the interaction picture they are replaced with 0 and the component a being acted on by the operator is evolved in Fourier space using a enc if you write a This means that lt CDATA dtheta_dz Lilthetal i theta theta theta theta damping dphi_dz L2 theta i phi phi phi phi damping no 1 gt the component theta will be evolved in Fourier space using the sum of the operators L1 and L2 and the evolution of phi will not include L2 theta Further writing something like 3 L1 theta does not have the intended effect theta will still only be evolved with L1 not 3 L1 We realise this is a potential source of error but it enables the equation syntax to remain uniform for all algorithms A future version of xmds will check the users equations more thoroughly When using explicit picture algorithms none of these problems exist If no errors are reported then what you get will be what your equations asked for Sometimes it is useful to calculate a function without having to calculate it separately for every transverse position This can be done by including the code in a lt functions gt element lt functions gt lt
45. gives modellers people doing computer modelling structure organisation and standardisation It provides a framework for describing the system they are trying to sim ulate be it in a physical mathematical scientific or even financial or economic setting It gives a way to keep the ideas behind a simulation well laid out and importantly docu mented for others to see and use And it gives people a common ground from which they can compare their numerical work something desperately lacking in an area at the interface between theory and experiment which already have a well ingrained culture of comparison and verification About this manual This manual has been split into five parts in an attempt to cover all of the material necessary to be able to use and master xmds but also to provide an entry point for novices and experts alike For the impatient it may be possible to learn jump right past the tutorial and learn what you need to get started from the worked examples in section 9 Part lis a very simple introduction to xmds and discusses how to build a simulation script from scratch how to modify existing code or a template to perform the simulation you want and how to perform stochastic simulations using the Message Passing Interface MPI Novice users may wish to start with Chapter 1 Starting from scratch or Chapter fl Using a template to help themselves get going with xmds Advanced users of xmds may 111 lv skip th
46. gt 11 22 3 5 2 vectors cross_propagation required vectors string variableName string variableName lt vectors gt Contains array of strings Subelement of Path to tag lt simulation gt lt sequence gt lt integrate gt lt cross propagation gt lt vectors gt Description The names of the cross propagating vectors as a list of strings Example lt simulation gt lt sequence gt lt integrate gt lt cross_propagation gt vectors main vc lt vectors gt lt cross_propagation gt lt integrate gt lt sequence gt lt simulation gt 11 22 4 breakpoint optional lt breakpoint gt xmds tags lt breakpoint gt Contains Subelement of Path to tag lt simulation gt lt sequence gt lt breakpoint gt Description Saves some vectors to a binary XSIL file when this element is reached in the simulation This can be used for generating an XSIL file that can be used as an input for another simulation or to check the behaviour of the simulation part way through This feature is described further in Chapter 2 Section Example 214 LANGUAGE REFERENCE lt simulation gt lt sequence gt lt breakpoint gt uL ds CARS gt lt breakpoint gt lt sequence gt lt simulation gt 11 22 4 1 filename breakpoint optional lt filename gt string lt filename gt Contains string Subelement of Path to tag lt simulation gt lt sequence gt lt brea
47. gt lt CDATA const double energy 4 const double vel 0 0 const double hwhm 1 0 ies lt globals gt lt field gt 135 136 WORKED EXAMPLE NLSE XMDS lt name gt main lt name gt dimensions t lt dimensions gt lt lattice gt 100 lt lattice gt domains 5 5 lt domains gt lt samples gt 1 lt samples gt lt vector gt lt name gt main lt name gt lt type gt complex lt type gt lt components gt phi lt components gt lt fourier_space gt no lt fourier_space gt lt CDATA const double wO hwhm sqrt 2 log 2 const double amp sqrt energy w0 sqrt M_PI 2 phi pcomplex amp exp t t w0 w0 vel t gt lt vector gt lt vector gt lt name gt vc1 lt name gt lt type gt double lt type gt lt components gt damping lt components gt lt fourier_space gt no lt fourier_space gt lt CDATA damping 1 0 1 exp pow t t 4 4 10 gt lt vector gt lt field gt lt sequence gt lt integrate gt lt algorithm gt RK4IP lt algorithm gt lt interval gt 20 lt interval gt lt lattice gt 1000 lt lattice gt lt samples gt 50 lt samples gt lt k_operators gt lt constant gt yes lt constant gt lt operator_names gt L lt operator_names gt lt CDATA L rcomplex 0 ktx kt 2 1 gt lt k_operators gt lt vectors gt main vcl lt vectors gt lt CDATA dphi dz L phi i phi phi phi phi damping 9 1 THE
48. gt lt simulation gt 11 21 2 dimensions optional lt dimensions gt string variableName string variableName lt dimensions gt Contains array of strings Subelement of Path to tag lt simulation gt lt field gt lt dimensions gt Description A space delimited array of the names of dimensions in the field other than the propagation dimension This element is therefore a list of the transverse field dimensions Example lt simulation gt lt field gt lt dimensions gt t lt dimensions gt lt field gt lt simulation gt 11 21 3 lattice field optional lt lattice gt int int lt lattice gt required if lt dimensions gt assignment present Contains array of integers Subelement of Path to tag lt simulation gt lt field gt lt lattice gt Description A space delimited array of integers giving the number of points in the grid for each dimension listed in the lt dimensions gt tag Exmaple 192 LANGUAGE REFERENCE lt simulation gt lt field gt lt lattice gt 100 lt lattice gt lt field gt lt simulation gt 11 21 4 domains optional lt domains gt double double double double lt domains gt required if lt dimensions gt assignment present Contains array of ordered pairs of doubles Subelement of Path to tag lt simulation gt lt field gt lt domains gt Description This tag specifies the domain range of each dimension listed in
49. gt 50 lt lattice gt lt moments gt P_out N_out lt moments gt lt CDATA Prout Ps N_out N Tes lt sampling gt lt group gt lt output gt lt simulation gt This simulation solves for the propagation of an optical pulse through a field of atoms having a transition frequency tuned to that of the optical pulse centre frequency The atoms are modelled as two level atoms The propagation equations shown in Equation 10 13 are deceptively simple OE t 2 a 95 OPA _ py Ot ONG 2 EP 10 13 The reality of this problem is that there are three components two propagating in the main propagation dimension t and the other in the transverse dimension z The component E is the electric field amplitude P the polarisation state of the atoms and N the excitation state of the atoms 1 being all in the ground state and 1 being all in the excited state 10 4 TLA XMDS 161 Lastly g is a coupling constant between the electric field and the atoms The curious feature of this set of PDEs in fact the very reason why it is chosen it as an example is that there exists a soliton solution for the electric field as shown in Equation 10 14 2 t E t z sech 10 14 gto to which is time lagged with propagation at the rate g a 10 15 The result for the electric field in the above simulation is shown in Figure 10 4 in which the soliton solution is evident
50. h F IN e FA a 11 ay e O7 ehe sti ab 12 x x az a3 h 13 ac a bsa b3 2a TENS e U a3 h x ki ac 15 a h F N 2 F7 a 16 a ea ehe sa a 17 z a a4 az h 18 ay a bia bia ba 39 19 ag e Capac ad 20 aa h F N x F aa 1 a4 hL 2 k 21 Ad el ad 22 x 7 as E a4 h 6 4 NUMERICAL METHODS FOR PARTIAL DIFFERENTIAL EQUATIONS 109 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 a a bsia bs a b5 39 bs Ag lt e as h 2 k Ae ae h F N 10 FT aJ eO ehe aL Ae x 22 ag as h a a 061a bea boa bo 4aq bs sa em e l a6 hL 2 ke f aj a h FIN 2 F aj amp as eU a6 h x k1 a x 22 a ag h a a b Aa b 2a b 32 b Ag b sa b a a e l ar h 2 ke1 aj aj h F N GET aj 0 aj eO he a k aj x 22 ag az h ay a bg 9 bsa bg 7a dos e as hC 29k 1 ab ay h FIN 25 F ay E Uas hL 20 k a x x ag ag h a a 091a bo 6a bora bg gay ss e U s h x ki ac a h F N 2 FT ac 8 T eO a9 hE 20 k ae x x ajo ag h 110 N
51. in the split step semi implicit method above This method was derived by Rob Ballagh s BEC group at the University of Otago It is described in detail in the PhD thesis of B M Caradoc Davies 7 We will present only an optimised recipe for implementing this algorithm 1 Calculate if required 2 ax a 3 a e Pe a 4 aj a 5 ax hN x ax 6 ax ehe akg T a a sax 8 39 x ih 9 ax ar jay 10 ax hN x ax ll a a iax 12 ax ar lax 13 ax hN x ax 14 a a lax 15 2 2 ih 16 ax a axK 17 ag ei 5k aK 18 ax hN x ax 19 a ezh a k1 a 20 a a ax 102 NUMERICAL MODELLING THEORY Similarly to the semi implicit algorithm in the interaction picture the first and the last steps are easily performed with the field in Fourier space since the linear operators reduce to ki multipliers And again note that the dependence of the linear operators Ci z9 ki on the propagation dimension ought to be weak in comparison to the size of the time step This algorithm has a similar memory overhead as its Explicit picture counterpart but it no longer needs the memory for the p derivatives vector nor does it have to do the work to calculate this vector The interaction picture method relies on the derivative Operator being independent of the propagation dimension otherwise it cannot function as a fourth order algorithm 6 4 11 The Fourth Fifth Order adaptive Runge Kutta Method in the
52. in xmds and people do in practice however one has to jump through some hoops that we don t really want to bother with here So to imitate Dirichlet boundary conditions we ll not let the solution evolve outside the domain of the transverse dimension z All this means is that we make that particular domain rather larger than was necessary and we make sure we don t evolve it in time for too long We also start with a very narrow initial condition and not have the diffusion coefficient too high so as to inhibit the diffusion a bit Please note that this is an example simulation to get people used to using the syntax of xmds and not necessarily to pedantically solve certain physical problems we merely use the physical situations here to illustrate xmds not the other way around An important solution of the diffusion equation is a Gaussian of the form 1 zr Faa where zo is the location of the maximum and the standard deviation c t increases with time as Tif 1 6 c t V2kt 1 7 This solution is also a handy initial condition and this is the analytical form of the initial condition we will be giving to xmds 22 STARTING FROM SCRATCH 1 2 2 Starting off the simulation code Now with the boundary conditions and the initial condition specified analytically we are now in a position to start writing some xmds code As per usual we start with the lt xml version 1 0 gt and lt simulation gt tags Then add the simu
53. k bo xmidS o v kf awk R mom ee ame amp 2 Q a xe AS ma dE e Re You 152 103 65bresumds 2 S Mok a aO O 28 o a a da BB Q 155 10 4 Ta ASA 158 ia SS SSS ee ee ee s ee eee 162 a Saas eee ma we e eS 167 IV Reference Manual 175 177 as ena a ye a a a a a 177 11 2 name imation ee cuu Su gum XE he eee awe aa 177 11 3 prop dim simulation ota eee eee ee eee aa ge d 178 LIA error C cile ee es Se ee en ee ed code 178 Gh Ae a Be Oe ee AAN ee 179 11 6 stochastic ss 179 11 7 MPI Method 180 See RN MEM a a RS NE as RETE M RES 180 Ll 9 seed vueae eoe Rer PR ese eC MAN de b G ue ere Ee 181 AA A 181 becaria dci a ee s 63k S ase 2548 e ee ee eoe 182 LLIZBmaty output 24 sa mes bo Rd 3e dela be BEES Xox doe 182 TT LICUIT IE 183 11 14use doublel e 183 ee te RM NER Sea S a NEN ESES SENE P NE S 184 LI 6threadsl soon ea m Ye ERE ue ie SORS e Red xe e Eq xo dede o 185 FEMENINA 185 DEE A A a A 186 a A ae ee ek ee a a A es Uni da de 186 ee MR 187 IEA Torpe TTL 187 11 20 2name 1818 2 Ee Bb ae ee sas s 188 11 20 3type Gb lle uy a p y Vene xal ice NRA eee BEE es 189 CONTENTS xi 11 20 4 default value 2 2 a 189 Ne MENG og dee s a xam sk ii Blok a XE kk a SE A 190 1121 1 name field 2 6 2 oe oe RE sus RE s n q Q a a ed Q s G W w Qua 190 11212dimensions s m a e a EE Rode Y SUQ 191 11 21 3 lattice Held 42444
54. libraries and library search paths required for using threaded fftw3 simulations 2 8 PREFERENCES 49 e User defined options at configuration USER_LIB if xmds has been installed in the user s home directory then this flag needs to be specified so that xmds can find the xmds specific libraries when building a simulation USER_INCLUDE if xmds has been installed in the user s home directory then this flag needs to be specified so that xmds can find the xmds specific header files when building a simulation 2 8 4 Examples of changing options Using gcc for g When configuring xmds before installation the configuration script often sets the default C C compiler for xmds to be g This is the GNU C compiler Sometimes this is not desirable and so one may wish to use the GNU C compiler gcc instead To do this one needs to change two things the XMDS_CC setting and the XMDS_LIBS setting In the xmds prefs file one then sets XMDS_CC to gcc and appends 1stdc to the list of flags already given for XMDS_LIBS at installation The addition of 1stdc is so that gcc can make use of the C extensions to gcc so that it can actually compile the simulation which is in some sense a mix of C and C Using icc for g An alternative C compiler to g is the Intel C Compiler icc To use this compiler instead of g or gcc assuming of course that you have the compiler installed on your system just set XMDS_CC to icc and pre
55. lt dimensions gt transverse dims gt 1 attice gt lt lattice gt SIL mGa OE pointer or SAC dam gt lt domains gt lt domains gt lt domain of each dimension gt lt s lt v lt lt fi lt amples gt lt samples gt ASES amp Mes cp Ono Sr ector gt lt name gt main lt name gt lt type gt complex lt type gt SUS cera typeof vector lt components gt lt components gt names of components gt X fourier space fourier space lt defined in k space gt lt CDATA ILES vector eld The sequence of integrations to perform sequence i ntegrate gt algorithm lt algorithm gt lt RK4EX RK4IP ARK45EX ARK45IP SIEX SI lt iterations gt lt iterations gt lt default 3 for SI algs gt lt interval gt lt interval gt AO Temasimedinmi Mo lattice lt lattice gt lt I mos polares alin madim Gam gt lt samples gt lt samples gt lt no pts in output moment group gt lt k_operators gt lt constant gt yes lt constant gt res mo lt operator_names gt lt operator_names gt lt CDATA 11 gt lt k_operators gt lt vectors gt lt vectors gt lt l esto memes gt lt CDATA 54 USING A TEMPLATE 1 gt lt integrate gt lt sequence gt lt hie output to generate gt lt output format binary precision single gt lt group gt lt
56. moment groups are being evaluated one being integrated over the t dimension and the other integrated over both transverse dimensions y 10 0 z z a Time integrated intensity b Total energy vs z FIGURE 10 1 Results for ndparamp xmds 10 2 kubo xmds STan u e r S O al 0 tee lt pele Kio Oscilleror slmilarlon 2 gt lt simulation gt lt name gt kubo lt name gt lt prop_dim gt t lt prop_dim gt lt error_check gt yes lt error_check gt lt stochastic gt yes lt stochastic gt lt paths gt 1 lt paths gt lt use_mpi gt no lt use_mpi gt lt seed gt 1 2 lt seed gt lt noises gt 1 lt noises gt lt field gt lt samples gt 1 lt samples gt lt vector gt lt name gt main lt name gt 10 10 2 KUBO XMDS 153 lt type gt complex lt type gt lt components gt z lt components gt lt CDATA ALES 11 lt vector gt lt field gt lt sequence gt lt integrate gt lt algorithm gt SIEX lt algorithm gt lt interval gt 10 lt interval gt lt lattice gt 1000 lt lattice gt lt samples gt 100 lt samples gt lt iterations gt 3 lt iterations gt lt CDATA dz dt i z n 1 11 lt integrate gt lt sequence gt lt output gt lt group gt lt sampling gt lt moments gt realz lt moments gt lt CDATA realz z gt lt sampling gt lt group gt lt output gt lt simulation gt The kubo oscillator is described in Equation 10
57. of stochastic calculus are Ito calculus and Stratonovich calculus Ito calculus is better suited to mathematical proofs but Stratonovich calculus has the major advantage that derivatives follow the normal chain rule This typically simplifies the numerical technique for solving a Stratonovich equation with a high order method Both types of equations arise naturally from modelling physical systems Fortunately it is possible to transform a set of equations in Ito calculus to Stratonovich calculus and vice versa If a set of Ito equations is given by zw salto lale 6 10 then the corresponding set of Stratonovich equations is O gat r e z o belgas abes cole 610 Numerical integration of stochastic equations is more computationally intensive than deter ministic equations of equivalent size and complexity 6 3 Numerical Methods for Differential Equations We now briefly describe several numerical methods for integrating deterministic and stochas tic differential equations 84 NUMERICAL MODELLING THEORY 6 3 1 The Euler and Inverse Euler Methods Consider the ordinary differential equation a aut flt a 6 12 where f is some functional of t and or a Given that a current value of a say a a tn is known then the simplest method to propagate a forward in time is to use the current slope flin an to calculate the new value a 1 An 1 An h f t anl 6 13 where A is the time step so that
58. on the xmds web site This is intentional as we set w 0 which should give the same answer however this simulation is more general than the later one discussed in Chapter The main difference careful readers will have noticed is that the average solution of z does not decay away to zero This is because we ve used a perturbation constant of 0 1 as opposed to Z used in the later example 4 1 2 Making the simulation hard Now s the time to up the ante In some real life situations lots and lots of paths are required so that decent statistics are generated For instance in modelling the evolution of stock prices one can use a stochastic differential equation however for banks etc to be really 4 2 WITH MPI 73 0 9 7 0 8 7 0 7 7 0 6 7 0 5 7 0 4 7 FIGURE 4 2 Scilab generated plot of the evolution of the oscillator frequency z over time perturbed by noise sure of the results they need to run of the order of millions of paths This can take a very long time especially if run on a single processor Therefore in such situations being able to automatically parallelise ones simulation and run the different paths over many CPUs would be great and this is something that xmds does really well as you shall hopefully see Let s run the simulation over a few more paths so we ll just glibly throw a some zeroes at the lt paths gt tag and set it to 1024000 and run the simulation again kubo tutorial Beg
59. operating systems including the Cygwin environment on Windows though an experienced programmer might easily port it to another platform xmds requires a very small amount of disk space and the RAM required is dependent algorithm but is usually little more than the size needed to store a certain number of copies of the fields to be integrated A C compiler is required as are also 126 FUNCTIONALITY BEGIN Repeat for Half and Full step Integration Loop for N Paths Initialise Segment 1 Sample Segment 2 Post Processing Add Path pi Compute Means and Standard Deviations Corhpute Half Full Step Error mie Output END FIGURE 8 1 xmds a functional diagram the FFTW Fastest Fourier Transform in the West libraries for C Further for stochastic problems xmds can produce output code which uses MPI routines to parallelise the problem for running on multiple CPUs or computer clusters An MPI compiler is required to take advantage of this XMDS can also use an OpenMP threads to use multiple processors on single jobs if an OpenMP enabled compiler is available 8 1 1 Installation Installation instructions and files are available at http www xmds org downloads First if parallel processing is desired a working MPI system and or an OpenMP enabled compiler should be installed Second the FFTW libraries must be installed if they are not already present These libraries may be f
60. particular dimension is in Fourier space ie if the corresponding lt fourier_space gt assignment was yes then it is referenced by its name as specified in the lt field gt element but prefixed by a k Finally be wary of unintentional integer divisions for example const double c 1 2 as the compiler will perform an integer division on the 1 2 before converting to a double resulting in c 0 0 It ought to be written as const double c 1 0 2 if c 0 5 was the intended result Worked example nlse xmds XML is very much like HTML so readers with any experience in writing HTML will find XML fairly straight forward However even if a user has not had any experience in HTML they will pick up the XML with the help of some example simulation files A comprehensive guide to general XML syntax is available in 9 For now a good introduction will be to step through the nlse xmds example script then we can work through some more advanced examples to elaborate on some of the options available These examples scripts can be found in the examples directory along with many others They may also be downloaded from the web site http www xmds org examples html Sod VOTTEN il Q ee lt gt Example simulta aon Non ner Schroedinger Equation z lt simulation gt lt name gt nlse lt name gt lt prop_dim gt z lt prop_dim gt lt error_check gt yes lt error_check gt lt stochastic gt no lt stochastic gt lt globals
61. phi 5 111 lt moment_group gt lt moment_group gt lt moments gt py ic lt moments gt lt integrate_dimension gt yes yes no lt integrate_dimension gt lt CDATA py phi 1 ic phi 5 Jes moment group moment group lt moments gt ppy ichi lt moments gt lt integrate_dimension gt no yes yes lt integrate_dimension gt 10 6 HIGHDIM_VECTOR_VERSION XMDS 171 lt CDATA ppy phi 1 ichi phi 5 11 moment group vectors main vc lt vectors gt lt CDATA const complex densi phi 2 phi 1 phi 4 phi 3 const complex dens3 const double gdensi const double gdens3 phi 5 re phi 5 re phi 5 im phi 5 im phi 6 re phi 6 re phi 6 im phi 6 im dphi_dt 1 L2p phi 1 i Vr 2 gamloss ii gamloss132 2 gamloss12 2 dx dy dz phi 1 i gameff gamlossi2 densi phi 1 i gam13 gamloss132 dens3 phi 1 i phi 2 phi 3 ixchippy ippy dphi_dt 2 L2n phil 2 i Vr 2 gamloss11 gamloss132 2 gamloss12 2 dx dy dz phi 2 ixgameff gamlossi2 densi xphi 2 i gam13 gamloss132 dens3 phi 2 i phi 1 phi 4 dphi_dt 3 L4p phi 3 i Vr 3 gamloss31 gamloss132 2 gamloss32 2 dx dy dz phi 3 i gam33 gamloss32 densS3 phi 3 i 0 5b phi 1 phi 1 i gami3 gamlossi32 densi phi 3 dphi_dt 4 L4n phi 4 i Vr 3 gamloss31 gamlossi32 2 gamloss32 2 dx dy dz phi 4 ixgam33 gamloss32 dens3x ph
62. problem Example lt simulation gt lt sequence gt lt f inter lt moment_group gt lt moments gt joe is great lt moments gt lt integrate_dimension gt no yes lt integrate_dimension gt lt CDATA joe psi psi is phi phi great theta gt lt moment_group gt lt filter gt lt sequence gt lt simulation gt 11 22 3 2 functions filter optional lt functions gt CDATA lt functions gt Contains CDATA Subelement of Path to tag lt simulation gt lt sequence gt lt filter gt lt functions gt Description This is the best place to define any functions that do not depend on the transverse dimensions Example 11 22 SEQUENCE 211 lt simulation gt lt sequence gt lt filter gt lt functions gt lt CDATA Some C code Wile lt functions gt lt filter gt lt sequence gt lt simulation gt 11 22 3 3 vectors filter required lt vectors gt string variableName string variableName lt vectors gt Contains array of strings Subelement of Path to tag lt simulation gt lt filter gt lt vectors gt Description The names of the vectors xmds should apply the filter to Given as a space separated list of strings Defaults to main Example lt simulation gt lt filter gt lt vectors gt main vcl lt vectors gt lt filter gt lt simulation gt 11 22 3 4 fourier space filter optional lt fourier_space gt boolean
63. sample 200 points and hence set lt samples gt to 200 The code for this looks like 12 STARTING FROM SCRATCH lt interval gt 10 lt interval gt lt lattice gt 10000 lt lattice gt lt samples gt 200 lt samples gt Having told xmds some of the parameters it must use to perform the integration we haven t yet told it how to actually carry out the integration By this I mean that we have to describe in terms of C language code the differential equation that xmds is to use to evolve the solution forward see Equation 1 3 We do this by using a CDATA block and writing the equation in a form understandable to both us and xmds which isn t technically speaking C code lt CDATA dx dt sipmat y x Ql oh xu Y IUE dz_dt x y bxz gt Notice that this looks similar to the analytical form of Equation in that we have described the derivatives with respect to time t of the fields z t y t and z t With that we have completed the lt integrate gt element and the lt sequence gt section The simulation script is now lt xml version il ONE lt simulation gt lt name gt lorenz lt name gt lt author gt Paul Cochrane lt author gt lt description gt Lorenz attractor example simulation Adapted from the example in Numerical methods for physics by Alejandro L Garcia page 78 1st ed lt description gt lt prop_dim gt t lt prop_dim gt lt globals gt lt CDATA const dou
64. short flag unique xmds chooses the next character in the variable name in this case the letter c and making the short form flag c The long form remains the same i e nconst Of course if the letter c has also been used then xmds loops through all of the other letters in the variable name until it finds a match If no match is found xmds reports an error To check the short and long forms of the flags the simulation expects the user can use the h or help options with the simulation For instance with the nlse simulation one can enter the following command to get a listing of the arguments the xmds simulation expects at its command line nlse help Example lt simulation gt lt argv gt lt arg gt lt name gt nconst lt name gt lt type gt int lt type gt lt default_value gt 2 lt default_value gt 11 20 ARGV 189 lt arg gt lt argv gt lt simulation gt 11 20 3 type arg required lt type gt char lt type gt xmds 1 2 Contains string Subelement of Path to tag simulation lt argv gt arg type Description The type of the command line argument e g int double char At present this mechanism cannot handle complex inputs Example simulation lt argv gt lt arg gt lt name gt nconst lt name gt lt type gt int lt type gt lt default_value gt 2 lt default_value gt lt arg gt lt argv gt lt simulation gt 11 20 4 defa
65. should use to do so We do this using the lt moments gt tag and a CDATA block For this simulation things are quite simple since we just want the amplitude of the variables as they evolve Just to make this really obvious we ll define the output moments to be new variables called x0ut yOut and zOut which are just equal to the variables x y and z but it makes it more obvious in the code what information we are grabbing out of xmds It is possible to define as many moments as you wish but you must define at least one The XML code we use is lt sampling gt moments x0ut yOut zOut lt moments gt lt CDATA LODEL yOut y zOut z 11 gt lt sampling gt giving an output section which looks like 1 1 A BASIC SIMULATION 15 lt output gt lt sampling gt lt moments gt x0ut y0ut zOut lt CDATA xDut X yDut y zOut Z gt lt sampling gt lt output gt lt moments gt and an overall simulation which is finally lt xml version 1 0 gt lt simulation gt lt name gt lorenz lt name gt lt author gt Paul Cochrane lt author gt lt description gt Lorenz attractor example simulation Numerical methods for physics page 78 1st ed in lt description gt lt prop_dim gt t lt prop_dim gt lt globals gt lt CDATA const double sigma 10 0 const double b 8 0 3 0 const double r 28 0 const double xo 1 0 initial const double yo 1 0 i
66. size of the output xmds since xmds 1 2 has had the ability to generate binary output files which are inherently smaller and can have better precision than ascii data files but also deals away with the redundancy introduced in the way that the ascii data is stored In xmds 1 2 binary output was controlled by the lt binary_output gt and lt use_double gt tags This syntax is deprecated as of xmds 1 3 in favour of passing attributes to the lt output gt element and this is the syntax we ll be discussing here Users who are still using xmds 1 2 are advised either to upgrade to a more recent version of xmds to make use of the better syntax If you still wish to use xmds 1 2 then the syntax for using binary output is described in Chapter 11 2 5 1 The format attribute As mentioned above the output format of an xmds simulation is now controlled by the format attribute of the lt output gt element The syntax is as follows lt output format ascii binary gt d other xmds tags gt lt output gt where by saying ascii binary we mean that the format option is either the string ascii or binary the double quotes are necessary and that the default option is ascii For those who have used xmds in the past you may remember that all of the data is output into the xsil file Binary output doesn t stop the generation of the xsil file but 38 EXTRA AND ADVANCED FEATURES merely uses a feature of the XSIL format that en
67. sqrt e2 2 M PI r2 r2 149 150 MORE EXAMPLES 1 gt lt globals gt lt field gt lt name gt main lt name gt lt dimensions gt y t lt dimensions gt lt lattice gt 100 100 lt lattice gt domains 10 10 10 10 lt domains gt lt samples gt 1 1 lt samples gt lt vector gt lt name gt main lt name gt lt type gt complex lt type gt lt components gt ff1 ff2 sh lt components gt lt fourier_space gt no no lt fourier_space gt lt CDATA ffl pcomplex ampirexp pow y yc1 r1 2 2 pow t tc1 r1 2 2 vy1x y ff2 pcomplex amp2 exp pow y yc2 r2 2 2 pow t tc2 r2 2 2 vy2 y sh rcomplex 0 0 iles lt vector gt lt vector gt lt name gt vci lt name gt lt type gt double lt type gt lt components gt damping lt components gt lt fourier_space gt no no lt fourier_space gt lt CDATA damping 1 0 i1 exp pow y y t t 8 8 10 gt lt vector gt lt field gt lt sequence gt lt integrate gt lt algorithm gt RK4IP lt algorithm gt lt interval gt 10 lt interval gt lt lattice gt 500 lt lattice gt lt samples gt 50 50 lt samples gt lt k_operators gt lt constant gt yes lt constant gt lt operator_names gt Lapi Lap2 lt operator_names gt lt CDATA Lapi Lap2 i ky ky kt kt 1 i ky ky kt kt 2 1 10 1 NDPARAMP XMDS 151 1 gt lt k_operators gt lt vectors gt main vci lt vectors gt l
68. that we used 16 computers instead of 1 and now we ve got a wall time improvement of about a factor of 4 As I mentioned above some of this will be due to communication overhead i e too much of the computer s CPU is used in just sending the data back and forth across the network and not enough time is actually spent just doing the calculation However we have reduced the amount of time necessary to come to a solution so we have actually done what we set out to do A new option for the MPI routines within xmds is to implement load balancing This is where the work of a parallel simulation is allocated to each of the CPUs running part of the job on the basis of the load of the CPU CPUs with less load on average get more work to do and CPUs with more load on average have less work to do The adaptive scheduler that performs this operation can be switched on by the use of the code lt MPI_Method gt Scheduling lt MPI_Method gt or disabled which may be useful if your system doesn t allow threads for example by the code lt MPI_Method gt Uniform lt MPI_Method gt After you ve finished your simulation and generated the results all that is left to do is to clean up after having run a parallel simulation You do this by running the wipe program Z wipe v lamhosts And that s it We ve managed to decrease the amount of wall time necessary to calculate many paths of a stochastic simulation and all it took was for us to change one li
69. the lt dimensions gt tag Example lt simulation gt lt field gt lt domains gt 5 5 lt domains gt lt field gt lt simulation gt 11 21 5 samples field required lt samples gt int int lt samples gt Contains array of integers specifically either 1 or 0 Subelement of Path to tag lt simulation gt lt field gt lt samples gt Description Tells xmds which moment groups to sample This tag should contain a space separated list of 1 s and 0 s A 1 meaning sample the moment group in question and a 0 meaning do not Example 11 21 FIELD 193 lt simulation gt lt field gt lt samples gt 1 lt samples gt lt field gt lt simulation gt 11 21 6 vector required lt vector gt xmds tags lt vector gt Contains filename lt type gt lt components gt lt fourier_space gt lt vectors gt CDATA Subelement of Path to tag lt simulation gt lt field gt lt vector gt Description A container for tags describing a vector of the field Example lt simulation gt lt field gt lt vector gt SSS samele EYEE 2 lt vector gt lt field gt lt simulation gt 11 21 6 1 name vector required lt name gt string lt name gt Contains string Subelement of Path to tag lt simulation gt lt field gt lt vector gt lt name gt Description Name of the vector in question At least one vector must be present and at least
70. the following relevant define statements placed at the head of the file field main defines define _main_ndims 1 define main latticeO 100 define main xminO 5 000000e 00 define main dxO 5 000000e 00 5 000000e 00 double 100 define main dkO 2 M_PI 5 000000e 00 5 000000e 00 define dt main dxO define dkt main dkO vector main defines define main main ncomponents 1 define phi main main main main pointer 0 vector sgi coterms defines 121 define _main_sgi_coterms_ncomponents 1 define _sgi_coterm_L_phi _main_sgi_coterms _main_sgi_coterms_pointer 0 This does not make for compact code but it was desired that the generated code remain readable by the user so that they may inspect it to see exactly what xmds is doing in regards to solving their problem It also has the added advantage that should a user wish to do something rare and obscure they only need to write a short script to generate a base C code listing which they may then alter to suit their problem 122 DEVELOPMENT AND PROGRAM STRUCTURE Functionality This chapter describes the full functionality of xmds but in practice most of this detail is not required to integrate common types of ODEs and PDEs One of the particularly useful features of xmds is that both the high level description of the problem and the generated program can efficiently handle systems with widely varying complexity Skip to the worked examples to see how v
71. this algorithm to second order there is no guarantee of convergence for stochastic equations 6 3 5 The Adaptive Fourth Order Runge Kutta Method For problems that exhibit very diverse behavior of the solution say with slow change for large parts of the propagation and rapid variations in others it is very useful to employ an algorithm that can adjust its stepsize This allows the simulation to speed through smooth uninteresting countryside in a few great strides but to tiptoe with many small steps through treacherous terrain The adaptive stepsize method implemented in xmds in an embedded fourth fifth order Runge Kutta method ARK45 It takes advantage of the fact that it is possible to combine certain six function evaluations that can result in a fifth order Runge Kutta method in another way to give a result that is accurate to fourth order The ARK45 algorithm calculates the forth and fifth order solution and uses the difference between them as an estimate of the current discretisation error The structure is ki A f tn an ko h f tn azh an b21k l k h f t azh an b31K b32ko ky h f t a4h an bak bagke ba3k3 k h f tn ash an bsiki b52k3 b53k3 b54k4 ke h f t agh an bork be2K2 be3k3 beaks besks Anti An ciki coke c3k3 c4k4 csks coke O A9 at a ctki ko cjk3 cika cks c kg O h 6 26 n Here a
72. this example is adapted from that in Garcia 1 We shall use in this example an initial pulse which is a cosine modulated Gaussian A x t 0 cos k z 2p exp A l 3 2 20 where k 2z is the wave number of the pulse which has a wavelength A zo is the initial peak position of the pulse and c is the initial pulse width standard deviation This initial 51 52 USING A TEMPLATE condition has an analytical solution which we give here merely for interest value A z t cos k x vt xo exp lr me ol j 3 3 cos k x zo vt exp f Ec wot ce 3 4 which is still a cosine modulated Gaussian just displaced by an amount vt The boundary conditions we shall use are periodic and since we get this for free when using xmds there isn t anything special we need to do here other than to keep in mind that the boundary conditions are periodic As with the diffusion equation example we ll make the x domain go from 0 5 to 0 5 we ll set the value of c to 0 1 and zo to 0 The wave speed we ll set to 1 Any other parameters will be mentioned as we go through and hack with the code 3 2 The template code Here is the code that we re going to hack with This is just an outline of many of the tags that we can use to perform a simulation using xmds Most of them you should know by now but one or two you may be unfamiliar with this is ok as we don t need them for what we want to do Comments
73. to for convergence of the method Defaults to 3 Example lt simulation gt lt sequence gt lt integrate gt lt iterations gt 5 lt iterations gt lt integrate gt lt sequence gt lt simulation gt 11 22 2 4 tolerance required lt tolerance gt double lt tolerance gt Contains double Subelement of Path to tag lt simulation gt lt sequence gt lt integrate gt lt tolerance gt Description If one of the adaptive step size methods is used this assignment controls the maximum relative error that is allowed per step on any of the grid points Example lt simulation gt lt sequence gt lt integrate gt lt tolerance gt ie 10 lt tolerance gt lt integrate gt lt sequence gt lt simulation gt 11 22 SEQUENCE 201 11 22 2 5 max iterations optional max iterations int max iterations Contains integer Subelement of Path to tag simulation sequence integrate max iterations Description Applies to the ARK45 algorithms only If the behaviour of the solution is unknown and might result in a very small step size the maximum number of iterations steps and discarded steps can be limited with this optional tag This is particularly useful in the debugging stage or on systems that impose time limits on their running programs Example lt simulation gt lt sequence gt lt integrate gt lt max_iterations gt 5000000 lt max_iteratio
74. to use complex here when one is using a Fourier transform technique to integrate the equations even if initially the variables are double Example lt simulation gt lt field gt lt vector gt lt type gt complex lt type gt lt vector gt lt field gt lt simulation gt 11 21 6 4 components required lt components gt string variableName string variableName lt components gt Contains array of strings Subelement of 196 LANGUAGE REFERENCE Path to tag lt simulation gt lt field gt lt vector gt lt components gt Description Names the components of the vector If a large array of components are required it is possible to declare that here by putting the size of the array in parentheses immediately after the name Subsequent references to these components must have the index afterwards in parentheses Example lt simulation gt lt field gt lt vector gt lt components gt phi JoelsGreat 256 lt components gt lt vector gt lt field gt lt simulation gt 11 21 6 5 fourier_space vector optional lt fourier_space gt bool bool lt fourier_space gt required if lt filename gt not present Contains array of booleans Subelement of Path to tag lt simulation gt lt field gt lt vector gt lt fourier_space gt Description Tells xmds whether the dimension specified is initialised in Fourier space This is a space separated list of yes no valu
75. used to build simulations e XMDS_INCLUDES the include flags used for building simulations e FFTW_LIBS the libraries specific to fftw e FFTW MPI LIBS the libraries specific to fftw but for the MPI compiler e FFTW3_LIBS the libraries specific to version 3 of fftw e MPICC the MPI C compiler e MPICCFLAGS the compiler flags to use with mpicc e MPILIBS the library flags to pass to mpicc e USER LIB the location of the file xmds library functions only necessary if xmds is compiled for user use e g L home cochrane bin e USER_INCLUDE the location of the xmds header files only necessary if xmds is com piled for user use e g I home cochrane bin Note when fftw is installed in usr local then with fftw path must be specified so that the configure script can find the libraries and headers 128 FUNCTIONALITY One can also set various options on the command line as arguments to the configure script The options available are e with user to install xmds into the user s bin directory within their home directory This will also point the USER_LIB and USER_INCLUDE variables to the correct places so that xmds will find its required headers and libraries so that it can build simulations e enable mpi to check for MPI in the configuration and to enable the use of MPI for building simulations e enable threads to check for the ability of the fftw libraries to use threads and enable their use in simulations e
76. when you have two variables to be entered at the command line that start with the letter k What xmds does to solve this problem is if a variable already has a short option taken e g if we had already defined another variable in the lt argv gt list called say kruntsch then the next character is used for the short option which would be the letter a for kappa Of course if this letter is taken then xmds searches for a single character representation of kappa throughout the variable name until it finds one that isn t used If xmds doesn t find a short option that isn t used then it throws an error Assuming that everything has worked ok and the assignments to the short options have worked properly how can one find out what the short option is if it has changed Well you can simply ask the simulation for help Just run the simulation with either h or help and it will print out the usage of the simulation and a list of the option names their data 46 EXTRA AND ADVANCED FEATURES type and default value For instance asking the diffusion simulation for help we get the following output diffusion help Usage diffusion k lt double gt Details Option Type Default value k kappa double 0 1 So we call diffusion with k and the simulation is expecting a double precision number after the k flag Also we are told that either k or kappa are possible options but we already knew that anyway and that kappa is
77. where the rate of change of the solution varies a lot or when the appropriate step size is unknown As with the fourth order Runge Kutta there is no guarantee of convergence for stochastic equations 6 4 Numerical Methods for Partial Differential Equa tions The methods in the previous section all deal with ordinary differential equations quite well but partial differential equations are considerably more complex since the derivative depends also upon the values for a at other points in the transverse space In principle once the transverse spaces have been divided into a finite lattice sets of partial differential equations are just larger sets of ordinary differential equations but there are many subtleties involved with achieving this discretisation The critical issue is how the evolution due to the transverse derivatives is performed This issue divides algorithms into two classes explicit picture and interaction picture We begin here by explaining just how the transverse derivatives may be calculated and then outline how the Semi Implicit and the Runge Kutta algorithms can be applied in each of these cases 6 4 1 Evaluating Transverse Derivatives The transverse dimensions must be discretised to form a lattice Lattices with more points enable more detail but also take up more memory and more computational resource The transverse derivatives must be evaluated at each lattice point and there must be some method for doing so tha
78. with the information you received as to the offer to distribute corre sponding source code This alternative is allowed only for noncommercial distribution and only if you received the program in object code or executable form with such an offer in accord with Subsection b above The source code for a work means the preferred form of the work for making modifications to it For an executable work complete source code means all the source code for all modules it contains plus any associated interface definition files plus the scripts used to control compilation and installation of the executable However as a special exception the source code distributed need not include anything that is normally distributed in either source or binary form with the major components compiler kernel and so on of the operating system on which the executable runs unless that component itself accompanies the executable If distribution of executable or object code is made by offering access to copy from a designated place then offering equivalent access to copy the source code from the same place counts as distribution of the source code even though third parties are not compelled to copy the source along with the object code 4 You may not copy modify sublicense or distribute the Program except as expressly 234 GNU GENERAL PUBLIC LICENSE provided under this License Any attempt otherwise to copy modify sublicense or distribute the Program is vo
79. 0e 00 0 000000e 00 0 500000e 00 1 000000e 00 1 000000e 00 0 000000e 00 1 000000e 00 0 000000e 00 0 000000e 00 1 000000e 00 1 000000e 00 0 000000e 00 0 000000e 00 1 000000e 00 0 500000e 00 0 000000e 00 0 000000e 00 1 000000e 00 0 000000e 00 0 000000e 00 0 000000e 00 1 000000e 00 0 500000e 00 0 000000e 00 0 000000e 00 1 000000e 00 1 000000e 00 0 000000e 00 0 000000e 00 lt Stream gt lt Array gt lt XSIL gt The xsil2graphics utility program Applications xmds xmds devel examples xsil2graphics help xsil2graphics utility supplied with xmds version 1 6 5 Usage xsil2graphics options infile Options infile required the input xsil file h help optional display this information m matlab optional produce matlab output default 8 8cilab optional produce scilab output a mathematica optional produce mathematica 5 x ouput e mathematica optional produce mathematica ouput g gnuplot optional produce gnuplot ouput r R optional produce R ouput o outfile optional alternate output file name v verbose optional verbose output For further help please see http www xmds org Most users will find this utility essential it writes the output as dat ascii files and writes a separate file which is executed from within a common plotting program in order to load in the data Currently supported packages include matlab http www mathworks com scilab http scilabsoft inria fr Mathemati
80. 2 k1 F aa aa N x aa P 98 NUMERICAL MODELLING THEORY 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 9T 58 59 60 61 62 63 64 65 x x ay ayo h a at bis ia bua b11 79 b11 884 b11 99 011 084 p F 2 k1 F ae ae N x ae P E 20 x aiz a11 h a a dy21aq bio 68 bi2 za bio sa biooa b12 1084 b12 1196 p F L 2 ki F az a N x ay p 20 x aig ar2 h ay a bisa bisa b13 7a dis ga b13 9a bis 104 b13 1192 bis 1291 p F L 2 k1 F aj a N x ay p 20 x aja a13 h a a b14 1aa dia ai b14 72 014 585 dia oac bia 1084 b14 11ae b14 12a f 0b14 139g p F L 2 k1 F an an N x an P E a at by51a bisa b15 7a bi5 gas b15 92 b15 108a4 bis 1192 dis 1207 b15139y sr bisa a 1 b16 6 b15 6 a b161 b15 1016 6 b156 92 016 6 b156 a 016 7 b15 7b16 6 015 6 a5 b16 8 b16 5b16 6 b15 6 a4 b16 9 b15 9b16 6 b15 6 9 b16 10 b15 10b16 6 b15 6 aa b16 11 b15 11016 6 b15 6 a b16 12 b15 12b16 6 b015 6 a f b16 13 515 13016 6 15 6 ag 016 14 b15 14616 6 b15 6 an x a aa a14 h p F L 2 k F ail a N x ai p 1 z aie a15 h p F L 2 k1 F aj aj N x aj p 6 4 NUMERICAL MET
81. 2 4 o due oe ee a we ea 11 23 2 2 post propagation List of Figures 1 1 Three dimensional plot in Matlab of the trajectories of a Lorenz attractor Parameters used were o 10 b 8 3 r 28 with initial conditions of Lo 1 0 yp 1 0 and zo 20 0 18 1 2 Three dimensional plot in Scilab of the trajectories of a Lorenz attractor Parameters used were o 10 b 8 3 r 28 with initial conditions of zo 1 0 yo 1 0 and zo 20 0 4 eorr a 20 1 3 Three dimensional plot in Matlab of the diffusion of Gaussian pulse according to the diffusion equation Parameters used were amp 0 1 6 0 1 x9 0 30 1 4 Three dimensional plot in Scilab of the diffusion of Gaussian pulse according to the diffusion equation Parameters used were amp 0 1 0 0 1 x9 O 31 3 1 Three dimensional plot in Matlab of a cosine modulated Gaussian pulse ac cording to the advection equation Parameters used were v 1 0 0 1 zo 0 k c o Notice that with the periodic boundary conditions that the pulse moves off one side of the figure and re enters from the opposite side 3 2 Three dimensional plot in Scilab of a cosine modulated Gaussian pulse ac cording to the advection equation Parameters used were v 1 c 0 1 zo 0 k c o Notice that with the periodic boundary conditions that the pulse moves off one side of the figure and re enters from the opposite side
82. 66 RS 6404 ar rr ERR RU S Eo ss 191 11 21 4domains 192 11 21 5samples Hield s ros bee a ea eae eae 192 MeL 21 0 VECTOR ris Hae da sa de AB ee mode SR Sy Rs Hs He Sh dos 193 11 21 6 1 name VectoE lt e o oe eee Ae CRE GEER EE ne 193 11 21 6 2 filename vector cenas EROR eee A s 194 Ma IN 195 11 21 6 4 components x s e e so doge ge da d ue a he A 195 11 21 6 5 fourier space vector 2 hk e x eee 196 11 21 6 6 vectors vector J doe deo Sure eue foe r a 196 a a E sap ees et 197 pibes ba ee wy od we ee a e 197 rrr 198 TL TCI TP 198 11 22 2 2 interval A 199 11 22 2 3 iterations lt i sos po sun xo op 3G S Y Q w a Q ae ie 200 11 22 2 4 tolerance II 200 11 22 25 max ilterations ooa da a a de 201 11 22 2 6 min E steple x9 xx plu das AA Ro 201 11 22 2 T smallmemory Le LES e eee eee a 202 TAI ora a e rs EER BSG AES 202 11 222 9halt non EE 203 11 22 2 1 attice integrate sais a ds 204 11 22 2 1 samples intesrate ii o o bee ae eRe ess sd 205 bse ee Ge ee ae n ee Bey a 205 11 22 2 12 1vectors E operators caos ee Be Oe eae 206 FET Re ORO d s ad hee ded d aed 206 BY Deus i de Are qo 207 11 22 2 13moment group Integrate s su as aux esas spa 207 11 22 2 1 unctions Inteerate i a sde sta ne Du aD we munis 208 11 22 2 15vectors integrate lt erecta Gee a ee Sed 209 MP A EE s e Te una a NE 209 11 22 3 1 m
83. 70 an EI x F7 an 71 ar eC ma hc a A ar 72 a a bi5 19 bisa bisza bis ga bi5 ga bisioaq 015 186 bis 1997 b15 139y T b15 14Ap 73 aj 1 b166 b15 6 a b16 1 b15 1016 6 b15 6 9 b16 6 b15 6 a b16 7 b15 7b16 6 b15 6 a bi br6 8b16 6 b15 6 a4 b16 9 b15 9b16 6 b15 6 a b16 10 E b15 10016 6 b15 6 aa b131 b15 11016 6 b15 6 ae b16 12 D15 12b016 6 b15 6 a5 b16 13 d15 13016 6 015 6 9g 516514 b15 14016 5 015 5 an 74 z9 z9 a15 a14 h 75 a hF V ns F a 6 76 x x a16 ai5 h TT Bg PIN E 78 a a CA Cgap C98 C108q 01186 C128 f C13Ay C148 158j C168 The method requires 10 copies of the field like its explicit picture counterpart Note no rotations are required for steps 15 and 16 since they both occur at the end of the time step 108 NUMERICAL MODELLING THEORY 6 4 13 The Eighth Ninth Order adaptive Runge Kutta Method in the Interaction Picture The Eighth Ninth order Runge Kutta also requires a large number of Fourier transforms since there is no temporal symmetry in the algorithm However the improvement in con vergence may offset the extra time required particularly in stochastic problems 1 Calculate if required 2 Se Bah 3 a a 4 a x a 5 aa h F N 22 95 EN bd eee aa 7 x9 z9 a aj h 8 a a 051a ee e U a h x k1 Ab 10 a
84. CDATA double f sin 2 0 z This is a function of the propagation dim LE lt functions gt Conversely it is often desirable to reference a variable in the equations that is a function of the spatial coordinates but is otherwise constant The variable damping is exactly one of these Rather than re calculating it at every time step and lattice point what we did here was to 144 WORKED EXAMPLE NLSE XMDS define an additional field vector vc1 in the lt field gt element that could store this variable calculate 1t once when it is initialised and then reference it in the main equations Therefore we must tell xmds which vectors we wish to access in the equations This is done with the lt vectors gt assignment The reason that we do not make all vectors available everywhere by default is one of efficiency For example if a particular vector was initialised in x space then it would have to be transformed to k space in order to use it in the lt k_operators gt code If we did not need to use it there then the Fourier transform would have been a waste of CPU time Further the damping variable used here would have to have been declared as type complex in order to have been able to be Fourier transformed to k space to be available in the lt k_operators gt element where it wasn t needed causing a waste of memory as well Another possible function required to define the equations of motion is some integral across the transverse di
85. ECTORS FROM FILE 39 the lt filename gt tag within the lt vector gt element within the lt field gt element The syntax for this is lt field gt lt vector gt lt filename format ascii binary xsil gt lt enter the file name here gt lt filename gt lt vector gt lt field gt where ascii is the default option when the format attribute is not specified As of xmds 1 3 xmds has the ability to load binary as well as ascii data Which xmds should expect is given by the format attribute of the lt filename gt tag within the lt vector gt element Using binary input however doesn t significantly change how the data should be organised prior to loading into an xmds simulation If MPI is enabled xmds will only load into memory the appropriate part of the input file irrespective of the file format 2 6 1 Intialisation from an XSIL file As of xmds 1 5 3 xmds can initialise a vector from a moment group of an XSIL file produced by a lt breakpoint gt tag see Section or an lt output gt tag in xmds If you are generating the XSIL file from an lt output gt tag then the output moment group must meet a certain format for xmds to be able to understand how to load the file correctly If the file is generated from an lt breakpoint gt tag then this is taken care of for you if the variables have the same names in the two simulations For XSIL files generated from output moment groups the format of the XSIL file mu
86. HODS FOR PARTIAL DIFFERENTIAL EQUATIONS 99 66 ar a 67 a a Ca Cgap C98 C108q Ci1de C128 f C13Ay C148 C158 i C168 68 a aru Cag C885 coa Cj98q cHa Conf C138g Cihan C158 C18 This method requires 11 extra copies of the field An additional copy is required for the adaptive time step compared to the non adaptive ARK9EX since the approximate error in the step must be calculated This is achieved using the field stored in a which holds the 8th order solution analogous to the a in the ARK45 method Once again we should note the adaptive step implemented in this method is stochastically safe making it the best performing solver for stochastic problems 6 4 8 Interaction Picture Methods Interaction picture methods can only be applied to equations of the general form a5 2 x1 F 2 ki Fl x A x a x x 6 37 Although it is possible to define interaction picture algorithms in which the linear opera tors may have the general matrix form here they are restricted to the diagonal form fi Non diagonal operators require eigen vector value decomposition each time step and at each lattice point which becomes computationally expensive These operators may still have coefficients that depend on the the propagation dimension and an the Fourier space coordinates k The evolution of the field a is carried out in both normal space and Fourier space in al
87. ILOS O 3 Sb o lt simulation gt 11 7 MPI_Method optional lt MPI_Method gt string lt MPI_Method gt Contains string Subelement of Path to tag lt simulation gt lt MPI Method gt Description Defaults to Scheduling For stochastic simulations using MPI it tells xmds which method to use for splitting the paths between different processors The default Scheduling option is usually optimal but requires the use of threads If threads are unavailable the Uniform option should be used Example lt simulation gt lt MPI_Method gt Uniform lt MPI_Method gt lt simulation gt 11 8 paths optional lt paths gt int lt paths gt required if lt stochastic gt is yes Contains integer Subelement of Path to tag simulation paths Description The number of stochastic paths to integrate Example lt simulation gt lt paths gt 3 lt paths gt lt simulation gt 11 9 SEED 181 11 9 seed optional lt seed gt int int lt seed gt required if lt stochastic gt is yes Contains array of two integers Subelement of Path to tag lt simulation gt lt seed gt Description The seed to the random number generator Internally these are the seeds passed to the C routine erand48 to generate two independent random number gen erators hence the use of two integers Example lt simulation gt lt seed gt 1 2 lt seed gt lt simulation gt 11 10 noises o
88. If set to yes the output executable will compare all output moments on a point by point basis between the full and half step runs It will then write out to screen the maximum value of the discrepancy and the half step results will be used for the final output The lt stochastic gt assignment is also optional the default being yes It tells xmds whether the simulation uses Gaussian noise terms or not If lt stochastic gt is set to yes then three further assignments lt paths gt lt seed gt and lt noises gt become compulsory The kubo xmds and fibre xmds examples in Chapter cover stochastic simulations in more 9 2 THE GLOBALS ELEMENT 139 detail The child elements within the lt simulation gt element that xmds then looks for are lt globals gt lt field gt lt sequence gt and lt output gt 9 2 The globals element lt globals gt lt CDATA const double energy 4 const double vel 0 0 const double hwhm 1 0 ie lt globals gt This element is optional and is used to define any numerical constants that are useful to have globally available to all sections of code The lt CDATA gt container tells the XML parser to copy its contents without attempting to parse them The contents const double energy 4 const double vel 0 0 const double hwhm 1 0 are copied directly into the output C code and therefore must have correct C syntax See Section 8 3 for more on C syntax
89. Interaction Picture In the Fourth Fifth Order adaptive Runge Kutta algorithm it seems we would require too many Fourier transforms for this to be efficient Following 8 the function can however primarily be evolved in Fourier space and transformed into normal space for calculation of the A x a only The use of the interaction picture does not allow a reduction of the computational effort to the same extent as in the RK4IP algorithm but this method can still be vastly superior over the adaptive explicit picture method as it allows larger step sizes for equations containing certain derivative operators In the following the vectors a a are supposed to be initially in Fourier space Note that for the calculation of the N these are transformed into normal space 1 Calculate if required 2 ax a 3 a ont a k E 4 ay a 5 a a 6 ax h FN 2 FL ax Tare g ah sk ak 8 a at caK 9 a a cjak 10 2 x ash 11 ax ar b3ak 12 ax h FIN 2 F ax 13 aj 1 bsi ci ar 03i c1a bar 6 4 NUMERICAL METHODS FOR PARTIAL DIFFERENTIAL EQUATIONS 103 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 x a az a3 h Note co c 0 eas a hc a9 ay ay ay h FIN x5 FT aj E as a h zk Az p ar 1 ba c1 ar b
90. Multiple path stochastic problems such as this may be parallised using MPI routines If an MPI compiler was specified in the configuration step of installation this will often be mpicc then all that needs to be done is to toggle the optional lt use_mpi gt assignment to yes xmds will then place the appropriate MPI calls in the output code and compile it with the MPI compiler The executable should then be run through the MPI execution handler probably mpirun with the number of processors option supplied For example if 16 processors are available then the final command for execution would be mpirun np 16 kubo Note that whole paths are assigned independently to the processors so there is no benefit in specifying more processors than there are paths in the simulation It is not necessary for the number of processors to be a factor of the number of paths some processors will simply do one more paths than others Also note that for both of the semi implicit algorithms SIEX and SIIP an iterations assignment may be used within the integrate element to specify the number of iterations to use in the method refer Sections and 6 4 9 This assignment is optional and will default to three when absent The reason that the kubo oscillator is used as an example is that it has an analytic solution as shown in Equation 10 10 Figure 10 2 shows the results for a single trajectory and an averaged trajectory illustrating the expected behavi
91. SIMULATION ELEMENT 137 gt lt integrate gt lt sequence gt lt output gt lt filename gt nlse xsil lt filename gt lt group gt lt sampling gt lt fourier_space gt no lt fourier_space gt lt lattice gt 50 lt lattice gt lt moments gt pow_dens lt moments gt lt CDATA pow_dens phixphi 111 lt sampling gt lt group gt lt output gt lt simulation gt The above XML file describes to xmds how to write a program that solves the Nonlinear Schr dinger Equation in one dimension Equation 9 1 A A ol 9 1 ot 1397 Here is the single component complex field phi in the input file is the propagation dimension z 7 is the local time coordinate t and T is a damping field damping applied to absorb scattered radiation at the domain boundaries Figure shows what the output of this simulation should look like At the head of the file lt a mver sion il 0 m lt Example simulation Non Linear Schroedinger Equation gt Here the first line defines the file as an XML file xmds requires correct XML syntax but otherwise does not need to validate the file with a DTD Document Type Definition The second line is simply a comment line All comments begin with a lt and end with a Any characters are allowed in between save for the comment closing sequence gt The remainder of the file is covered in detail in Sections 9 1 9 8 while Section B 2 co
92. UMERICAL MODELLING THEORY 48 ag a bioiaa bio sa 197a biosa bio pa 49 ag AN 50 ag h F N 2 E lab 51 Ad eU a0 h x k1 i ad 52 x x a11 ajo h 53 ae a 041184 b11 59 dura bisa birga O11 108a 54 a e Iar hL 20 k a 55 a h FIN 25 F ae amp 56 ae eU an h xxki B ae 57 aq z9 a12 au h 98 ay a 05184 bio sa 01578 bio 884 bio oac 0151984 bio ia 59 ar e ei hE a9 M ar 60 ag h FIN 2 F las 61 ar e mes af 62 z9 z9 a13 az2 h 63 ag a bisa bisa b13 7a bis ga b13 92 b13 10 a b13 1192 bis 1291 64 ag e Uars hL 2 k ay 65 ay WF IN x FA a E 66 ay etre q 67 z9 x aa zs aj3 h 68 a at 514184 bisa 01478 bia say 014 99 014 10 4 b14 118e 0141297 biasa 69 a e ars hL 20 k 7 a 70 h F N x F7 an 71 a ec eiae ski ar 72 a a bisia bisa 015 785 15 88 bis ga 0151084 bisa Dis 1291 b15139y RU bisa 6 4 NUMERICAL METHODS FOR PARTIAL DIFFERENTIAL EQUATIONS 111 73 aj 1 b166 b156 a b16 1 b15 1b16 6 015 6 2 b16 6 015 6 a 016 7 b15 7b16 6 b15 6 a5 b16 3 b16 5016 6 b15 6 a4 b16 9 515 9016 6 015 5 ac 016 10 b15 110b16 6 b15 6 24 b16 11 b15 11616 6 b15 6 ae b16 12 b15 12b16 6 b15 6 a f b16 13 015 113016 6 015 5 ag 16 14 b15 14016 5 015 5 an 74 22 x a
93. _dim gt lt simulation gt 11 4 error check optional lt error_check gt bool lt error_check gt Contains boolean Subelement of Path to tag lt simulation gt lt error_check gt 11 5 USE_MPI 179 Description Whether or not to run the simulation at half the time step as well as at the full time step and give the difference between the results Defaults to yes Example lt simulation gt lt error_check gt yes lt error_check gt lt simulation gt 11 5 use mpi optional use mpi bool use mpi Contains boolean Subelement of Path to tag simulation use mpi Description Whether or not to use MPI routines for parallel processing of the simulation This writes very different code depending on whether the simulation is stochastic or non stochastic Only certain problems on certain systems can be efficiently parallelised when the equations are non stochastic Defaults to no Example lt simulation gt lt use_mpi gt yes lt use_mpi gt lt simulation gt 11 6 stochastic optional lt stochastic gt bool lt stochastic gt Contains boolean Subelement of Path to tag lt simulation gt lt stochastic gt Description Defaults to no Tells xmds whether or not the simulation uses Gaussian noise terms If this tag is set to yes then the lt paths gt lt seeds gt and lt noises gt tags become compulsory Example 180 LANGUAGE REFERENCE lt simulation gt ISO MS S
94. _space gt no lt fourier_space gt lt CDATA phi 0 es lt vector gt lt field gt lt sequence gt lt integrate gt lt algorithm gt SIIP lt algorithm gt lt interval gt 2 5 lt interval gt lt lattice gt 5000 lt lattice gt lt samples gt 50 lt samples gt lt k_operators gt lt constant gt yes lt constant gt lt operator_names gt L lt operator_names gt lt CDATA L i kx kx 11 gt lt k_operators gt lt iterations gt 3 lt iterations gt lt CDATA dphi_dt L phi ggamma phi beta sqrt 2 complex n_1 n_2 gt lt integrate gt lt sequence gt lt output gt lt group gt lt sampling gt lt fourier_space gt yes lt fourier_space gt lt lattice gt 50 lt lattice gt lt moments gt pow_dens lt moments gt lt CDATA pow_dens conj phi phi 135 lt sampling gt lt group gt lt output gt lt simulation gt This simulation solves Equation 10 11 in which a one dimensional damped field is subject to a complex noise This is a stochastic PDE OW y 0 uL amp z t i o z D 10 11 cg W Ola ies 10 11 Again the reason for using this as an example of a stochastic PDE is that it has an analytic solution as shown in Equation 10 12 Figure displays the results of this 158 MORE EXAMPLES simulation in Fourier space for a single trajectory and an averaged trajectory which appear as expected Idk t ao 4 where L is t
95. a double precision number of default value 0 1 And that s it At present xmds can accept int double float and char for command line arguments Complex numbers aren t yet implemented as of xmds 1 3 but may be added in a future version Now imagine that we wanted to run the diffusion simulation over a range of values starting from 0 1 to 1 0 To do this we could write a simple shell script as follows bin sh xgue a ati 0 1 0 2 OaS Oo4 0 5 0 6 0 7 0 8 0 9 1 0 do echo Running diffusion with kappa i echo i e diffusion kappa i diffusion kappa i mv diffusion xsil diffusion i xsil done Notice that we ve moved the output data file to a new filename since the file diffusion xsil will be produced each time the simulation it is run and hence our data it would get written over with new data each time the simulation is run had we not bothered to rename the xsil file Equivalently we could have used a Perl script to do the same thing for instance usr bin perl w use strict my i 0 1 while i lt 1 0 print Running diffusion with kappa i n printi are RU S TON kappa yal my Cargs diffusion kappa i system args mv diffusion xsil diffusion_ i xsil i i 0 1 Feel free now to extend and have a play with diffusion xmds For example change the 2 8 PREFERENCES 47 simulation to make it possible to vary sigma and even x0 and see how the output changes
96. a4 h ae a 05184 bs a b5 sac b5 au p F L O ki lal a N x ae P E x 22 ag as h a a bei 565a bsa b6 Ad besace p F L 1 k F az a N x ai p 20 x az ag h aj a b7 1aa b7 22 b 39 b7 484 b 59 b ea p F L 2 k1 F aj aj N x aj p z a ag a7 h a a bg 19 bg ga bs ra p F L 2 k1 F aj ay N x a p x a ag ag h ac a by a bo a bora bo 8a p J L 2 k1 F a a N X ac P x x aro e ag h 6 4 NUMERICAL METHODS FOR PARTIAL DIFFERENTIAL EQUATIONS 95 39 40 41 42 43 44 45 46 47 48 49 90 5l 52 53 54 99 96 57 58 59 60 61 62 ag a bip 19 bio 6a 510 78 bio sas b10 9ac p F L 2 k1 F ad ay N x aq p 2 x ay ayo h ae a buda bisa b11 79 b11 say b11 99 bis 1092 p F L 2 k1 Z ae ae N x ae P 2 x aig a11 h ay a bip aa 015 68 b12 7a bio 889 b12 9a 0121084 bio ia p F L 2 k1 F ay ar N x ap p x9 x aig ar2 h ag a b13 19 biasa bis a bis sas b13 9a bis 1oaa 013 186 dis 1997 p F L 25 k1 Z aj ay N x ay p 2 z9 aja a13 h a a b14 19 b14 69 b14 7a b14 8a b14 99 b14 10 4
97. ables binary files to be pointed to by the xsil file Therefore all of the important parameters of the simulation are still saved to the xsil file just the data is now saved to another file or files if you have more than one moment group containing just a binary string of data So when using binary output the following files will be produced a xsil file containing simulation parameters and pointing to the output data by default this will be called lt simulation name gt xsil and a binary data file for each moment group being called in general lt simulation name gt mg lt moment group number gt dat Running the atomlaser simulation with the format set to ascii we get an output xsil file of size 808 kB Now if we run the atomlaser simulation again except with the lt output gt tag set to lt output format binary gt then we get a xsil file of size 4 kB and a dat file called atomlasermg0 dat of size 336 kB giving a total of 340 kB which is 42 smaller than with just ascii output Bigger savings can be expected with longer simulations and or simulations using more dimensions 2 5 2 The precision attribute The default binary output is at double precision This is not always necessary for output of data especially if the data is to be displayed graphically and then interpreted further there the extra precision is not necessarily worthwhile Therefore there is also the precision attribute available in the lt output gt ele
98. ace into its respective Fourier domain that a partial derivative with respect to position just becomes i i e y 1 times the coordinate in Fourier space For our example such a mapping would be o gt tk 1 5 Az 1 5 Hence in Fourier space the second derivative on the right hand side of Equation 1 4 is just k Given that discrete Fourier transforms aren t that hard to do on a computer and in xmds we use the Fastest Fourier Transforms in the West http www fftw org so they re pretty fast using such a transformation improves the calculation somewhat 1 2 1 Specifying the problem We now need to specify the problem properly To do this we must specify an initial condi tion for the solution we wish to evolve and we must specify the boundary conditions of the domain over which we wish to solve this particular problem Boundary conditions are nec essary here since we have a transverse dimension i e x in this system recall in Section 1 1 we had no transverse dimensions only the propagation dimension of time In following Garcia 1 as we do here again we are trying to model the temperature diffusion of an initial temperature distribution in a one dimensional rod the ends of which are kept at a constant temperature of T 0 Unfortunately this implies Dirichlet boundary conditions and xmds only implements periodic boundary conditions This isn t strictly true as one can implement absorbing boundary conditions
99. ame gt lt type gt complex lt type gt lt components gt phi 6 lt components gt lt fourier_space gt no no no lt fourier_space gt lt vectors gt vcl lt vectors gt lt CDATA const double realfn mu 0 5xMxVr 1 Uoh11 hbar a Jrs J A NC 119222102259 phi j realfn gt 0 complex sqrt realfn 0 complex 0 0 else phi j complex 0 0 11 lt vector gt lt field gt 170 MORE EXAMPLES lt gt he siequence ot integrations to perform gt lt sequence gt lt integrate gt lt algorithm gt ARK89EX lt algorithm gt lt interval gt 1e 7 lt interval gt lt tolerance gt 1 0e 7 lt tolerance gt lt lattice gt 1000 lt lattice gt lt samples gt 10 10 1 lt samples gt lt k_operators gt lt constant gt yes lt constant gt lt operator_names gt L2p L2n L4p L4n lt operator_names gt lt CDATA L2p complex 0 hbar M 2 chix kx kx ky ky kz kz L2n complex 0 hbar M 2 chix kx kx ky ky kz kz L4p complex 0 hbar M 4 chix kx kx ky ky kz kz L4n complex 0 hbar M 4 chix kx kx ky ky kz kz gt lt k_operators gt lt moment_group gt lt moments gt chippy lt moments gt lt integrate_dimension gt yes yes yes lt integrate_dimension gt lt CDATA chippy phi 5 phi 5 Jes moment group moment group lt moments gt ippy ichippy lt moments gt lt integrate_dimension gt no no no lt integrate_dimension gt lt CDATA ippy phi 1 ichippy
100. arly with the RK4 method in both the explicit and interaction pictures the calcu lation of b is performed using the RK4 method and again using linear interpolation in the x dimension to create the mid lattice point values for a that are required Note that this use of simple linear interpolation causes a loss of order when implementing a higher order method such as the RK4 but generally the transverse lattice is sufficiently fine for this to be negligible However as mentioned earlier the semi implicit method in the interaction picture need only sweep through the memory space of main field vector once with each nonlinear step Therefore it would be desirable to use an algorithm that can calculate the cross vector si multaneously in the same sweep The semi implicit method in the interaction picture does 112 NUMERICAL MODELLING THEORY this by employing the standard semi implicit algorithm in the propagation direction simul taneously with another semi implicit algorithm in the transverse propagation dimension In pseudo code it looks like this hL 2 k 1 l a e a 2 Calculate if required 3 39 z9 h 4 For each point in space do a a a b br b c For N 1 iterations do i a a 3hN xa b ii b bz Ax H x a b d a 2 ar hN x a b e b b Az H x a b f b x Az br Ax H x a b 5 19 x h hL x k 6 a e a This routine assumes that the cross vec
101. asy to write down in xmds This is especially true since the derivatives are with respect to time for if we had spatial derivatives we would have to use Fourier transform techniques and mappings to simplify the calculation we ll cover this stuff in Section 1 2 so don t worry that we re not discussing it here which would complicate our script a bit more and we re trying to keep things simple here The Lorenz model can be used in the study of many interesting phenomena however it is possibly best known as a model of global weather 1 For a system to be chaotic it must be extremely sensitive to initial conditions such that any small perterbation of the initial 6 STARTING FROM SCRATCH conditions will cause wildly divergent evolution and this is something we will hopefully see here so let s continue There can be many global parameters in an xmds simulation although one in particular is special This is the propagation direction specified by the prop dim tag which is specified as part of the global functionality of the xmds simulation and appears within the simulation tags normally after the description In the problem we are trying to solve our field is evolving in time given by the variable t in the above equations and so the field is said to be propagating in t so we use time as the propagation dimension Therefore we add the line lt prop_dim gt t lt prop_dim gt to our xmds script just after the lt descript
102. ata contained in nlse xsil into matlab 229 230 LOADXSIL M UTILITY SCRIPT C 5 Authors Written by Paul Cochrane C 6 Bugs No known bugs However the loadxsil script does not work in Matlab version 4 0 or below it can only be used with Matlab version 5 0 and above Users with Matlab 4 0 can use the xsil2graphics utility as a means to import data into Matlab C 7 See also xmds 1 xsil2graphics 1 http www xmds org C 8 Copyright Copyright C 2003 2004 Code contributed by Paul Cochrane This program is free software you can redistribute it and or modify it under the terms of the GNU General Public License as published by the Free Software Foundation either version 2 of the License or at your option any later version This program is distributed in the hope that it will be useful but WITHOUT ANY WARRANTY without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE See the GNU General Public License for more details You should have received a copy of the GNU General Public License along with this program if not write to the Free Software Foundation Inc 59 Temple Place Suite 330 Boston MA 02111 1307 USA Gnu General Public License GNU GENERAL PUBLIC LICENSE Version 2 June 1991 Copyright C 1989 1991 Free Software Foundation Inc 51 Franklin St Fifth Floor Boston MA 02110 1301 USA Everyone is permitted to copy and distribute verbatim copies of this license document but
103. ate gt lt sequence gt lt simulation gt 11 22 SEQUENCE 209 11 22 2 15 vectors integrate required lt vectors gt string variableName string variableName lt vectors gt Contains array of strings Subelement of Path to tag lt simulation gt lt sequence gt lt integrate gt lt vectors gt Description The vectors xmds needs to access in the equations Given as a space sepa rated list of strings Defaults to main Example lt simulation gt lt sequence gt lt integrate gt lt vectors gt main vcl lt vectors gt lt integrate gt lt sequence gt lt simulation gt 11 22 3 filter optional lt filter gt xmds tags lt filter gt Contains fourier space COMA Subelement of Path to tag lt simulation gt lt sequence gt lt filter gt Description Container element for the tags describing how the field is to be filtered if at all Example lt simulation gt lt sequence gt lt filter gt lt l xmds tags gt lt filter gt lt sequence gt lt simulation gt 210 LANGUAGE REFERENCE 11 22 3 1 moment group filter optional moment group xmds tags lt moment group Contains Integrate dimension CDATA Subelement of Path to tag lt simulation gt lt sequence gt lt filter gt lt moment_group gt Description Defines and calculates a number or vector integrated through zero or more of the transverse dimensions of the
104. ate gt element the CDATA section within the lt integrate gt element is lt CDATA dA dt L A gt Note that we didn t write down the equations as lt CDATA dA_dt v L A gt and define the operator as L rcomplex 0 0 kx This is because of the way xmds transplants the equations into the code and specifying the constants within the operators is the most general way xmds can do this Heed a warning here though xmds at present doesn t pick this kind of error up so one can write the equation as just mentioned and xmds will happily do its thing however your answers will be wrong So the advice here is to be really careful The lt sequence gt element now looks like this 58 USING A TEMPLATE The sequence of integrations to perform 2 lt sequence gt lt integrate gt algorithm RK4EX lt algorithm gt lt RK4EX RK4IP ARK45EX ARK45IP SIEX SIIP lt interval gt 1 lt interval gt Sos O rar in mem chaine a lt lattice gt 500 lt lattice gt SUL nos Homes sli menst n Claim gt lt samples gt 50 lt samples gt lt no pts in output moment group gt lt k_operators gt lt constant gt yes lt constant gt lt l yes m gt lt operator_names gt L lt operator_names gt lt CDATA L rcomplex 0 0 v kx 1 gt lt k_operators gt lt CDATA dA_dt L A 11 lt integrate gt lt sequence gt 3 3 5 The output to generate We re
105. ation 10 6 ninj 10 6 in MAk Ar 10 6 where m are the transverse dimensions in Fourier space and n the transverse dimensions in normal space Also note that 2 Tmax Tmin Within the main integration equations the variances must also reflect the integration step size as given by Equation 10 8 ij ru 10 8 HA Az pu nin y xmds uses the Box Mueller technique as shown in Equation 10 9 which generates a pair of Gaussian noises and from a pair of random numbers x and z gt that have a uniform distribution between zero and one El 169 2AIn 2 ee Dia lt y y y 0 lt y lt l yeR 10 9 Estimating the error between full and half step sizes now poses an interesting problem since both evolutions must use the same underlying noise which is a function of both space and time if the difference between these paths is to be meaningful The random number generator must be reset before each of these integrations but how is the noise in the half step case appropriated There are two methods The first is to use the same noise for both half steps as is used for the one whole step This is undesirable since it makes sense to use the half step integration results being the more accurate for the final output and therefore independent noises must be used for each step So the second solution is to do just this and use the average of the two noises when calculating the full step i
106. backward plan Keeping accumulated wisdom Finished fftw calculations Beginning full step integration Sampled field for moment group 1 at t Sampled field for moment group 1 at t 0 000000e 00 2 500000e 09 2 5 BINARY OUTPUT 37 Sampled field for moment group 1 at t 5 000000e 09 Sampled field for moment group 1 at t 7 500000e 09 Sampled field for moment group 1 at t 1 000000e 08 Beginning half step integration Sampled field for moment group 1 at t 0 000000e 00 Sampled field for moment group 1 at t 2 500000e 09 Sampled field for moment group 1 at t 5 000000e 09 Sampled field for moment group 1 at t 7 500000e 09 Sampled field for moment group 1 at t 1 000000e 08 maximum step error in moment group 1 was 3 802825e 11 Time elapsed for simulation is 10 seconds You will note if you have run the atomlaser simulation both with and without wisdom how quickly the simulation starts once some fftw wisdom is used Also the simulation tells you that it is using previously generated wisdom and that it is saving it for future use 2 5 Binary output When performing big simulations i e over many dimensions or when propagating for a large distance over the propagation dimension one is going to produce very large output files This can be a problem and the problem will be exacerbated by the fact that by default xmds outputs data in ascii format with a lot of redundancy As a way to reduce the
107. bia b14 1297 biasa p F L 2 k1 F a an N x an pP E a at bis bisa b15 7a bi5 gas bis 92 bis 1084 bis 1192 dis 1297 b15 139y ur bisa a 1 516 5 015 5 a D16 1 015 1016 5 015 5 aa 016 6 015 5 amp i D16 7 015 7016 5 015 5 a5 b16 8 016 5016 6 015 5 av b16 9 515 9016 6 015 5 ac 016 10 015 10016 6 015 aa b16 11 b15 11016 6 b15 6 ae b16 12 b15 12b16 6 b15 6 a5 61613 b15 13016 6 b15 6 ag b16 14 b15 14616 6 b15 6 an x z9 a15 aj4 h p F L 2 k F ail a N x ai P 96 NUMERICAL MODELLING THEORY 63 z9 x a16 a5 h 64 p F L x9 k1 Z aj 65 a N x aj p 66 a a Cia Cgap C99 C108qd Cie C128 f C138g C148j C15A C168 This method requires 10 extra copies of the field a to jj Many of the a b and c constants are zero in the Eight Ninth order Runge Kutta method are zero meaning not all fields produced at each midstep are required all the time this allows us to reduce the total number of fields required significantly The 15th and 16th steps were merged together as well to reduce the total fields needed in memory Stochastic problems when solved with this methods will achieve a significant performance boost many orders or magnitude greater than extra memory use 6 4 7 The Eighth Ninth Order adaptive Runge Kutta Method in the Explicit Picture
108. ble sigma 10 0 const double b 8 0 3 0 const double r 28 0 const double xo 1 0 initial conditions const double yo 1 0 initial conditions const double zo 20 0 initial conditions 11 lt globals gt lt field gt lt name gt main lt name gt lt samples gt 1 lt samples gt 1 1 A BASIC SIMULATION 13 lt vector gt lt name gt main lt name gt lt type gt double lt type gt lt components gt x y z lt components gt lt CDATA X XO y ye ZENONE gt lt vector gt lt field gt lt sequence gt lt integrate gt lt algorithm gt RK4EX lt algorithm gt lt interval gt 10 lt interval gt lt lattice gt 10000 lt lattice gt lt samples gt 200 lt samples gt lt CDATA dx dt sigma y x Ghy_chg wise gt y 28 dz_dt x y bxz gt lt integrate gt lt sequence gt lt gt yet more xmds codes to Soma lt simulation gt 1 1 6 The output element Ok we re almost there This is the last section that we need to worry about So just to recap we ve told xmds the general features and variables it should use to construct the simulation the field to integrate over and the way in which the integration should take place and most importantly the differential equation that xmds should use to evolve the solution That sounds like about it doesn t it Well no We haven t told xmds to output anything yet and it s a bit silly to spe
109. boot on n9 becO9 1 CPU Executing hboot on n10 bec10 1 CPU Executing hboot on n11 bec11 1 CPU Executing hboot on n12 bec12 1 CPU 76 STOCHASTIC SIMULATIONS AND MPI Executing hboot on n13 bec13 1 CPU Executing boot on n14 bec14 1 CPU Executing hboot on n15 bec15 1 CPU topology done Now MPI is set up we need change one line of code to make our simulation parallel Change the lt use_mpi gt tag to read yes and you re done Ok now run xmds over the script to rebuild the binary executable giving you something akin to the following xmds kubo_tutorial xmds Output file name defaulting to kubo_tutorial xsil compiling for MPI parallel execution mpicc pthread 03 ffast math funroll all loops fomit frame pointer lm lmpi llam 1stdc I home cochrane bin L home cochrane bin o kubo_tutorial kubo tutorial cc I home cochrane bin I usr local include lfftw threads lfftw L usr local lib kubo tutorial ready to execute We now run the simulation by using the command Z mpirun c 16 kubo tutorial The command line flag c to the mpirun program tells mpirun how many processes to generate over the machines we have and since we have 16 computers to choose from we set this to 16 here it could be higher if you want to run more than one process per CPU Running this simulation took almost exactly 6 minutes Well that wasn t a great speed improvement was it Especially given
110. ca http www wolfram com gnuplot and R If no output file name is specified then the last extension is removed from the input file name and an appropriate ending m sci or nb is appended To avoid confusion between data from different XSIL data containers as there may have been more than one container within the input file specified xsil2graphics appends all variable names with an integer specifying 227 228 THE XSIL2GRAPHICS UTILITY PROGRAM which XSIL data container in sequence they came from There is an example plotting program for Mathematica in the examples directory To use the output of xsil2graphics in matlab one merely needs to call the name of the xsil file without the xsil extension For example for the nlse simulation the output file is nlse xsil and the command one uses in matlab is gt gt nlse In scilab for the same simulation the command to use is gt exec nlse sci loadxsil m utility script C 1 Name loadxsil load simulation data into matlab C 2 Synopsis loadxsil xsil file gt C 3 Description Utility script bundled with xmds used to load simulation output data from the xsil data file into matlab where the results can be presented graphically To load data from the xsil file data_file xsil enter at the matlab command prompt gt gt loadxsil data_file xsil C 4 Examples At the matlab command prompt gt gt loadxsil nlse xsil loads the d
111. ce code And you must show them these terms so they know their rights We protect your rights with two steps 1 copyright the software and 2 offer you this license which gives you legal permission to copy distribute and or modify the software 231 232 GNU GENERAL PUBLIC LICENSE Also for each author s protection and ours we want to make certain that everyone understands that there is no warranty for this free software If the software is modified by someone else and passed on we want its recipients to know that what they have is not the original so that any problems introduced by others will not reflect on the original authors reputations Finally any free program is threatened constantly by software patents We wish to avoid the danger that redistributors of a free program will individually obtain patent licenses in effect making the program proprietary To prevent this we have made it clear that any patent must be licensed for everyone s free use or not licensed at all The precise terms and conditions for copying distribution and modification follow GNU GENERAL PUBLIC LICENSE TERMS AND CONDITIONS FOR COPYING DISTRIBUTION AND MODIFICATION 0 This License applies to any program or other work which contains a notice placed by the copyright holder saying it may be distributed under the terms of this General Public License The Program below refers to any such program or work and a work based on the Program means eit
112. changing it is not allowed Preamble The licenses for most software are designed to take away your freedom to share and change it By contrast the GNU General Public License is intended to guarantee your freedom to share and change free software to make sure the software is free for all its users This General Public License applies to most of the Free Software Foundation s software and to any other program whose authors commit to using it Some other Free Software Foundation software is covered by the GNU Library General Public License instead You can apply it to your programs t00 When we speak of free software we are referring to freedom not price Our General Public Licenses are designed to make sure that you have the freedom to distribute copies of free software and charge for this service if you wish that you receive source code or can get it if you want it that you can change the software or use pieces of it in new free programs and that you know you can do these things To protect your rights we need to make restrictions that forbid anyone to deny you these rights or to ask you to surrender the rights These restrictions translate to certain responsibilities for you if you distribute copies of the software or if you modify it For example if you distribute copies of such a program whether gratis or for a fee you must give the recipients all the rights that you have You must make sure that they too receive or can get the sour
113. chi F g delta const double biggamma g g delta 2 const double gami3 Uoh13 chi 10 5 HIGHDIM XMDS 163 const double gam33 Uoh33 chi const double gameff Uoh11 biggamma chi const double gamlossii loss11 2 chi const double gamlossi2 loss12 chi const double gamloss31 1l0ss31 2 chi const double gamloss32 10ss32 chi const double gamloss132 loss132 chi const double cnoise noise sqrt 2 0 1 gt lt globals gt lt argv gt lt arg gt lt name gt kjoek lt name gt lt type gt double lt type gt lt default_value gt 1 0e6 lt default_value gt lt arg gt lt arg gt lt name gt joekappamax lt name gt lt type gt double lt type gt lt default_value gt 1 0e2 lt default_value gt lt arg gt lt argv gt lt Field to be integrated over gt lt field gt lt dimensions gt x y z lt dimensions gt lt lattice gt 16 16 16 lt lattice gt lt domains gt 1 2e 4 1 2e 4 1 2e 4 1 2e 4 8 0e 3 8 0e 3 lt domains gt lt samples gt 1 1 1 lt samples gt lt vector gt lt name gt vci lt name gt lt type gt double lt type gt lt components gt vcore Vir V3r gVir gV3r lt components gt lt fourier_space gt no no no lt fourier_space gt lt CDATA vcore omegax omegax x xtomegay omegay y ytomegaz omegaz z z Vir 0 5 M vcore hbar chi gameff gam13 2 2 dx dy dz V3r M vcore hbar chi gam13 2 gam33 2 dx dy dz gVir 0 5 M vcore hbar chi gV3r M vcore hbar chi
114. component whereas in the interaction picture linear operators are dealt with separately to nonlinear operators Once a set of partial differential equations has been discretised to a finite lattice in the transverse dimensions it is formally equivalent to a large set of coupled ordinary differential equations and therefore it can be solved by higher order methods such as the semi implicit method of Section or the fourth order Runge Kutta RK4 method of Section 6 3 4 Partial derivatives of the fields can be determined by Crank Nicholson methods or spectral methods Thus explicit picture methods are capable of solving the fully fledged generalised PDE ac N x a x p x x 6 34 pi x FA EL 09 k F ai x 6 35 This is the PDE that is listed in Equation with the cross propagating vector b omitted Refer to Chapter 8 for a full explanation of the variables The important thing to note is that the total derivative of the vector a is expressed in terms of the general nonlinear functionals M and that this form allows for nonlinear partial derivatives that may have spatial dependence for example i s 6 4 3 The Semi Implicit method in the Explicit Picture Here we use the Fourier Transform or Spectral method to determine the transverse deriva tives The procedure after Section 6 3 3 is 90 NUMERICAL MODELLING THEORY 8 Calculate if required 2 2 h aj a For N 1 iterations do a
115. cript to switch back and forth between MPI and non MPI code is a single flag at the top Set the lt use mpi gt tag to read yes or no as desired How to run the resulting program will depend on the specific implementations of MPI that you are using We will show a typical case based on the LAM MPI implementation another major implementation is MPICH but details will vary from system to system 4 2 1 Example using LAM MPI One of the first things to do is to make a hosts file so that LAM can know what machines it needs to send jobs to The hosts file is just a list of hostnames of computers you have access to that can run MPI For example here is a sample hosts file called lamhosts cat lamhosts bec00 becO1 becO2 bec03 bec04 becO5 bec06 becO7 bec08 bec09 bec10 bec11 bec12 bec13 bec14 bec15 As you can see there are 16 computers available to run jobs on Next it s a good idea to do a recon to check that all of your hosts are working and they can accept MPI connections To do this run the following command and you should see similar output recon v lamhosts recon testing n0 bec00 recon testing ni bec01 recon testing n2 bec02 recon testing n3 bec03 recon testing n4 bec04 4 2 WITH MPI 75 recon testing n5 bec05 recon testing n6 bec06 recon testing n7 bec07 recon testing n8 bec08 recon testing n9 bec09 recon testing n10 bec10 recon test
116. cription Defaults to 2 This tag tells xmds which version of the fftw library to use Currently fftw3 does not support distributed memory fourier transforms using MPI but only shared memory transforms using threads Consequently fftw version 3 can only be used for simulations that are not MPI using deterministic simulations Stochas tic MPI simulations can take advantage of fftw3 as they do not need distributed memory fourier transforms To override the libraries that xmds links against when compiling a simulation that uses fftw3 override the FFTW3 LIBS library in your xmds prefs file Example lt simulation gt lt tt tweens ion S fftwi views ion gt lt simulation gt 11 19 globals optional lt globals gt C code lt globals gt Contains CDATA block with C code Subelement of Path to tag lt simulation gt lt globals gt Description Defines variables and constants that are globally available to all sections of the output C code Example 11 20 ARGV 187 lt simulation gt lt globals gt lt CDATA const double energy 4 const double vel 0 0 const double hwhm 1 0 gt lt globals gt lt simulation gt 11 20 argv optional lt argv gt xmds tags lt argv gt xmds 1 2 Contains Subelement of Path to tag lt simulation gt lt argv gt Description Overall tag containing the arguments to be supplied at the simulation com mand line Command line arguments can be very u
117. ctionality in the output from xmds simulations The code now is lt xml version 1 0 gt lt simulation gt lt name gt lorenz lt name gt lt author gt Paul Cochrane lt author gt lt description gt Lorenz attractor example simulation Adapted from the example in Numerical methods for physics by Alejandro L Garcia page 78 ist ed lt description gt 1 1 A BASIC SIMULATION 5 lt gt yet more xmds codes O come lt simulation gt Alternatively you can add such information to the script by putting comments in the source code The above code may then look something like this lt xml version 1 0 gt S empleo slmolerlome Women gt lt Adapted from the example in Numerical methods for physics gt SUS ly Ae jandro IL Camelia pa sal lt mes script ys Paul Cochirans gt lt simulation gt lt name gt lorenz lt name gt A yet more xmds code to come c lt simulation gt 1 1 2 General simulation options Now that the very basic preliminaries are out of the way we need to focus on the problem we are trying to solve So for the sake of argument let s try to solve the Lorenz equations dx m o u x 1 1 dy 4a TE Y a 1 2 dz p pe bz 1 3 where c r and b are positive constants and the variables have the initial conditions x t 0 zo y t 0 y z t 0 20 Even though we are trying to solve something as complex as a chaotic system this is very e
118. cture 101 6 4 11 The Fourth Fifth Order adaptive Runge Kutta Method in the Inter Soe RE A A a ae ae we A MAA 102 6 4 12 The Ninth Order Runge Kutta Method in the Interaction Picture 104 6 4 13 The Eighth Ninth Order adaptive Runge Kutta Method in the In teraction Picture ECKE Ee ES se e s 108 6 4 14 Adding a Cross Vector 111 6 5 Discretisati n and Sampling Errors amp 4 2 oo dus 64 444 da ut oe 112 III User Manual 115 7 Development and Program Structure 117 8 Functionality 123 8 1 Installing and Running xmds 4 4 44 Ro ox Re Rn a m nes 125 A du Q s s LIN nn ISTE E Kee Cadran 126 B2 UBAEB x nose s doge S ESR Gode hg ok ode e AA A d 128 8 1 3 Preferences sak wk obse ke pas Re eR a ee 129 8 2 Syntax SIA uos due aes A RR Ole See oe A aod 130 8 3 C codine within elements amp 2 5333 424 444 das X xe E v sg 132 x CONTENTS 135 9 1 The simulation element 137 9 2 The globals l ment uode e amp de ed Ue dox E A OS RUE de 139 9 3 The field element 139 9 4 The vector element 140 9 5 The sequence element ad Lilas pal o OE a A oe oS 141 9 6 Th integrate element o usines 44 5002 68 40 484 464 bo ok x 142 9 7 The filter element ll RR 144 9 8 The output element 22 oa m Ro no 4 ko RR GR RAE sexe 145 149 XCTI 149 10 2
119. d by the two goals of preserving the free status of all derivatives of our free software and of promoting the sharing and reuse of software generally NO WARRANTY 11 BECAUSE THE PROGRAM IS LICENSED FREE OF CHARGE THERE IS NO WARRANTY FOR THE PROGRAM TO THE EXTENT PERMITTED BY APPLICA BLE LAW EXCEPT WHEN OTHERWISE STATED IN WRITING THE COPYRIGHT HOLDERS AND OR OTHER PARTIES PROVIDE THE PROGRAM AS IS WITHOUT WARRANTY OF ANY KIND EITHER EXPRESSED OR IMPLIED INCLUDING BUT NOT LIMITED TO THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FIT NESS FOR A PARTICULAR PURPOSE THE ENTIRE RISK AS TO THE QUALITY AND PERFORMANCE OF THE PROGRAM IS WITH YOU SHOULD THE PROGRAM PROVE DEFECTIVE YOU ASSUME THE COST OF ALL NECESSARY SERVICING REPAIR OR CORRECTION 12 INNO EVENT UNLESS REQUIRED BY APPLICABLE LAW OR AGREED TO IN WRITING WILL ANY COPYRIGHT HOLDER OR ANY OTHER PARTY WHO MAY MODIFY AND OR REDISTRIBUTE THE PROGRAM AS PERMITTED ABOVE BE LIABLE TO YOU FOR DAMAGES INCLUDING ANY GENERAL SPECIAL IN CIDENTAL OR CONSEQUENTIAL DAMAGES ARISING OUT OF THE USE OR IN ABILITY TO USE THE PROGRAM INCLUDING BUT NOT LIMITED TO LOSS OF DATA OR DATA BEING RENDERED INACCURATE OR LOSSES SUSTAINED BY YOU OR THIRD PARTIES OR A FAILURE OF THE PROGRAM TO OPERATE WITH ANY OTHER PROGRAMS EVEN IF SUCH HOLDER OR OTHER PARTY HAS BEEN ADVISED OF THE POSSIBILITY OF SUCH DAMAGES END OF TERMS AND CONDITIONS 236 GNU GENERAL PUBLIC LICENSE Part VI Bibliogra
120. defaulting to diffusion sci Proccessing xsil data container 1 Writing data container 1 to file 1 2 7 1 Matlab Loading the data into Matlab or Octave gt gt diffusion Doing a whos Name Size Bytes Class error_temperature_1 50x51 20400 double array t_1 1x51 408 double array temperature_1 50x51 20400 double array x 1x50 400 double array Grand total is 5201 elements using 41608 bytes Plotting the data 30 STARTING FROM SCRATCH gt gt mesh t_1 x_1 temperature_1 gt gt be Orb gt gt ylabel x gt gt labe L u 7 you should see a figure similar to Figure e 2 15 1 NN ANS S IKKS 0 5 W DANOS SO OO SOS NES S 03 So 0 LOSSES gt gt SSSSSSSSSSS SSS SSSESSSESSSSSS pe SSSSS FIGURE 1 3 Three dimensional plot in Matlab of the diffusion of Gaussian pulse according to the diffusion equation Parameters used were amp 0 1 0 0 1 zo 0 1 2 7 2 Scilab Loading the data into Scilab gt exec diffusion sci gt temp_di zeros 50 51 gt t 1 zeros 1 51 gt temp_d2 zeros 50 51 x 1 zeros 1 50 temperature 1 zeros 50 51 error temperature 1 zeros 50 51 gt diffusioni fscanfMat diffusioni dat Error Info buffer is too small too many columns in your file 1 2 A MORE COMPLEX SIMULATION 31 gt temp_di diffusion1 1 temp d2 diffusioni 2
121. e Kutta Method in the Interaction Picture The ninth order Runge Kutta does require a large amount of Fourier transforms since there is no temporal Symmetry in the algorithm However the improvement in convergence may offset the extra time required particularly in stochastic problems 1 Calculate if required 2 39 o bah 3 ay a 4 a e a 5 aa h F N GP a ME NS chc aa 7 39 x ag a4 h 8 ay a 051a 1 a2 hL 2 k 9 ay e ab 10 a h FIN 2 F ay 6 4 NUMERICAL METHODS FOR PARTIAL DIFFERENTIAL EQUATIONS 105 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 NE eU a2 h x k1 ab x z9 az az h a a bz a 535a dcs e U 03 h x ki ae ac h FIN F ac a eU a3 hL x k1 a x 22 a4 az h aq a ba 124 ba gay ba 39 aic e laa h 2 ke ad aa h F N 2 F aa E ag el whe 2 ki x 22 ag ay4 h ae a bs 9 bs 2a bs sac bs 4Ag A e7 as hL 1 k ae a h F N 2 F ae Hem eO as he x91 TA x x ag as h a a bsa b6 2ap b6 3a b6 4aa basa a e as he x des a a h FIN 2 F aj amp eU a6 h x k1 a x a az ag h a a b a b7 22 b7 39 b7 484
122. e equations are in this general form they may be propagated in the first order dimension x using one of the algorithms detailed further on 6 1 1 Boundary Conditions In the case of partial differential equations the equations must apply either over all of the space spanned by the transverse dimensions or else over a confined region of this space Either way computational resources cannot solve over an infinite region of space and so boundaries of some sort must be imposed Therefore a description of what happens at these boundaries is essential for solving the problem There are three common types of boundary conditions 1 Dirichlet boundary conditions This is where the value of a component is fixed at a particular boundary pr Se 0 6 4 2 Neumann boundary conditions This is where the cross boundary gradient of a component is fixed at a particular boundary al la cu E 6 5 3 Periodic boundary conditions This is where the value of a component and all of its cross boundary derivatives is the same at each end of the dimension in question It is as if this dimension was originally circular except that it has been cut for the purpose of discretisation DN E pr ae ax vei os gx A n 0 12 69 If an infinite region of space is absolutely necessary then it can be mapped to a finite region of space using a transformation such as x tan e which maps an infinite range in x to a finite range 5
123. e of loose mode is that it allows one to break up a simulation into parts where each part requires a slightly different grid For example in the diffusion example in Chapter 1 Section the restriction was made that the simulation is not evolved for long enough such that the field becomes non zero at the edge of the grid With loose mode after running the simulation for some time on a small grid if the state of the field is sampled at the end of the simulation the simulation can be continued on a larger grid though still keeping the same step size in that dimension and ensuring that the grid points do overlap Also if one wishes to increase or reduce the number of points in a given dimension and keep the width constant initialise the state of the field with that dimension in k space as in this case the requirement that the step size in the k space dimension be the same is equivalent to the requirement that the width of that dimension in z space remain the same Hence the number of points in z space in that dimension can be increased or reduced Note that a binary XSIL file produced on any architecture can be used on any other architecture byte swapping is automatically done if the endianness of the machine running the simulation is different to the endianness of the XSIL file and XSIL files with the output in single precision can also be used In summary the syntax for initialisation of a vector from an XSIL file is field l
124. e of the Program is restricted in certain countries either by patents or by copyrighted interfaces the original copyright holder who places the Program under this License may add an explicit geographical distribution limitation excluding those countries so that distribution is permitted only in or among countries not thus excluded In such case this License incorporates the limitation as if written in the body of this License 9 The Free Software Foundation may publish revised and or new versions of the General Public License from time to time Such new versions will be similar in spirit to the present version but may differ in detail to address new problems or concerns 235 Each version is given a distinguishing version number If the Program specifies a version number of this License which applies to it and any later version you have the option of following the terms and conditions either of that version or of any later version published by the Free Software Foundation If the Program does not specify a version number of this License you may choose any version ever published by the Free Software Foundation 10 If you wish to incorporate parts of the Program into other free programs whose distribution conditions are different write to the author to ask for permission For software which is copyrighted by the Free Software Foundation write to the Free Software Foundation we sometimes make exceptions for this Our decision will be guide
125. e_prefs gt xmds 1 24 Contains boolean Subelement of Path to tag simulation use prefs gt Description Defaults to yes Tells xmds whether or not to use the user preferences file called xmds prefs located in either the user s xmds directory or within the directory local to the xmds simulation script By default xmds will use the preferences file if it exists if not then it uses the preferences set when xmds was built One can explicitly use the build preferences by setting use prefs to no The format of the preferences file is a sequence of key value pairs separated by an equals sign For instance if one wished to change the compiler used by xmds to compile simulations then one needs to change the XMDS_CC variable Hence to do set this to icc for example one would enter the following line into the xmds prefs file XMDS CC icc Example lt simulation gt lt use_prefs gt no lt use_prefs gt lt simulation gt 11 16 THREADS 185 11 16 threads optional lt threads gt int lt threads gt Contains integer Subelement of Path to tag lt simulation gt lt threads gt Description If set to a value of n 1 this tag tells xmds to create a simulation that uses the threaded version of the FFT library and to use n threads for calculating the FFT s Example lt simulation gt lt threads gt 3 lt threads gt lt simulation gt 11 17 use openmp optional lt use_openmp gt
126. eature of xmds is for you In versions of xmds before xmds 1 2 to be able to map a parameter space or run the program over many different values of a simple global variable you had to modify your script rerun xmds with its implied compliation step and then run the simulation for each value This put plainly sucked so we put in a way to pass arguments to the simulation binary executable enabling us to write a simple shell script or Perl or Python to run our simulation over many different values This removed the need to recompile the simulation again and again and generally speaking speeds things up and takes at least some of the pain out of doing things like mapping parameter spaces So how do we tell xmds to make the simulation accept command line arguments You do this with an lt argv gt tagset which you put somewhere before the lt globals gt element For those of you who have worked with C before and passing arguments to programs will notice that we ve used the argv name here for the list of arguments the program will accept in exactly the same way that C programming does by convention To set up this list we need to specify the arguments and the relevant properties of the arguments As such we need to tell xmds what the name of the argument is its data type and its default value for the instances when we don t want to specify the value on the command line As might be obvious here we have a nested structure of information
127. ement may have as many of the other sub elements as desired to perform the calculation and the child lt sequence gt s 1 1 A BASIC SIMULATION 11 can contain another element the lt cycles gt element which controls how many times a given lt sequence gt is repeated The lt cycles gt element is optional and defaults to one It is important to note that the order of segments specified within a sequence are significant and operations given will be performed in that order So to summarise most of the time you will just use one lt sequence gt and it will usually only contain just the one lt integrate gt section hence the code will look like lt sequence gt lt integrate gt lt More mes tags to come gt lt integrate gt lt sequence gt We shall try to discuss the other features and tags in more depth later on in more advanced tutorials 1 1 5 The integrate element Since the lt integrate gt element is quite complex and it does all of the hard work we ll spend some time discussing it We need to tell xmds the algorithm to use to integrate the field specified earlier To do this we use the lt algorithm gt tag This tag is optional and will default to SIEX for stocastic simulations and to RK4EX for non stochastic simulations however it is a very good idea to explicitly specify what algorithm your simulation is using if only to help yourself in six months time or a colleague who may end up readin
128. ences between Matlab or Octave and Scilab and this is why we discuss the two here XSIL files can also easily be translated into scripts suited for input into Mathematica gnuplot or R using the bundled software Before we can start using Matlab or Scilab we must convert the data contained in the xsil file into something that Matlab Octave or Scilab can understand To do this we use the utility program bundled with xmds called xsil2graphics To generate an input file for Matlab or Octave use either xsil2graphics lorenz xsil or xsil2graphics matlab lorenz xsil however the second example is redundant as a Matlab or Octave m file is the default output from xsil2graphics You should see in the current directory a file called lorenz m and a data file for the one moment group that we sampled for lorenzi dat For Scilab use xsil2graphics scilab lorenz xsil giving the files lorenz sci and lorenzi dat Ok now we re ready to fire up our relevant numerical processing and graphical environment and visualise the results 1 1 7 1 Matlab and Octave Start Matlab or Octave and once at the command prompt load the information contained in the data file by using the command gt gt lorenz doing a whos should give you something similar to this gt gt whos Name Size Bytes Class error_x0ut_1 1x201 1608 double array 18 STARTING FROM SCRATCH error_y0ut_1 1x201 1608 double array error_z0ut_1 1x201 1608 double array t_1 1x201 1608
129. ents gt assignment This isn t going to be defined in Fourier space so the lt fourier_space gt assignment is set to no The trickiest part here is now defining the equation in the CDATA block Recalling the initial condition in Equation 3 2 we therefore declare the variable A as A rcomplex cos k x x0 exp x x0 x x0 2 0 xsigma sigma 0 0 where we have used the rcomplex function with a zero imaginary argument to define the initially real quantity however of complex type The result of these definitions gives us the following code for the lt field gt element lt gt Eteld to be intesrated over gt gt gt lt field gt lt name gt main lt name gt lt dimensions gt x lt dimensions gt lt transverse dims gt lt lattice gt 50 lt lattice gt SL MO jes 20 Garcia claim gt lt domains gt 0 5 0 5 lt domains gt lt domain of each dim gt lt samples gt 1 lt samples gt lt io pample iler Poimp Ow Chimp gt lt vector gt lt name gt main lt name gt lt type gt complex lt type gt lt l datali type ofi vector gt lt components gt A lt components gt I names jor components c lt fourier_space gt no lt fourier_space gt lt def in k space gt lt CDATA A rcomplex cos k x x0 exp x x0 x x0 2 0 sigma sigma 3 3 RIPPING IT TO BITS 57 11 gt lt vector gt lt field gt 3 3 4 The sequence of integrations to perform N
130. ery different kinds of equations can be easily integrated xmds is designed to integrate the following general PDE ral N x a x p x b x gtx f de f x a x p bC E 3 8 1 p x FT EL CAT F a x 8 2 Z be H x a x b x 8 3 where the vector a represents an m component real or complex valued field though for some problems the field may only have one component The vector x is the real valued space in which a lies which is the propagation dimension 2 plus the transverse dimensions 27 if any The number of transverse dimensions that the field may have is limited to 64 by xmds but is more likely to be further limited by the complier and stack size used when the output code is compiled The vector p may have any number of components including zero and may be derived through the action of a matrix of linear operators on the main vector components The action of these linear operators is calculated with the main vector in Fourier space hence as explained in Section transverse partial derivatives reduce simply to multiplication by the transverse Fourier space dimensions The vector b is optional and is an n component real or complex valued field which propagates along the transverse dimension z For stochastic problems is a vector of independent real Gaussian noises The total derivatives for vectors a and b are expressed as the general functionals W and H respectively Note that this form allows fo
131. es A value of yes meaning the dimension is defined in Fourier space and no meaning the dimension is defined in x space Example lt simulation gt lt field gt lt vector gt lt fourier_space gt no lt fourier_space gt lt vector gt lt field gt lt simulation gt 11 21 6 6 vectors vector required lt vectors gt string variableName string variableName lt vectors gt Contains array of strings 11 22 SEQUENCE 197 Subelement of Path to tag lt simulation gt lt field gt lt vector gt lt vectors gt Description Tells xmds the names of the variables are to be referenced in a CDATA block Defaults to main Example lt simulation gt lt field gt lt vector gt lt vectors gt main vc lt vectors gt lt vector gt lt field gt lt simulation gt 11 22 sequence required lt sequence gt xmds tags lt sequence gt Contains lt integrate gt lt filter gt Subelement of lt simulation gt or Path to tag lt simulation gt lt sequence gt lt sequence gt Description Container tag for the sequence of integrations to perform May contain other sequences within itself If subsequences exist then they must contain a lt cycles gt assignment Example lt simulation gt lt sequence gt SL Suda tags gt lt sequence gt lt simulation gt 11 22 1 cycles optional lt cycles gt int lt cycles gt required in nested lt sequence g
132. ese chapters and may be more interested in how to use xmds with MPI for running parallel simulations which is discussed near the end of Part Part Tl is a reasonably general discussion of numerical techniques for solution of differen tial equations Some of the techniques discussed are used within xmds and the procedures used internally in xmds are outlined Part is an overview of the language and contains some of the details of working with and the workings of xmds This along with Part is the essence of the xmds 1 0 documentation In time this Part will be updated but at present only superficial changes have been made in its layout Part TV gives specific information about the keywords used in the xmds language their place within a script s structure examples of usage and what xmds expects as arguments within the respective tags This Part is particularly aimed at users familiar with xmds who are likely to want to know the syntax of a particular keyword or other details of a particular keyword Part V is an appendix and covers the XSIL output format the xsil2graphics utility program the loadxsil m script the GNU General Public License and a bibliography Tools used to build xmds These are the multifarious tools with which xmds its documentation both handmade and automatically generated and its web pages has been made e GNU development suite autoconf automake etc e General development tools make gcc g cvs e Edi
133. ess non finite numbers slower than finite numbers others process them faster Regardless they are unlikely to be useful results from an integration Note that by default xmds compiles the generated C code with optimizations that sacrifice strict compliance with floating point arithmetic standards e g IEEE 754 Many of these optimizations assume that all floating point numbers are finite Natu rally this may pose problems for a non finite number check In some cases the check for non finite numbers is optimized to oblivion Depending on other optimizations the time step may become NaN and the integration will never move forward and therefore never terminate An easy way to see if the check for non finite numbers survived optimization is to run the following supposing xmds compiled the simulation binary from simulation xmds strings simulation grep halt_non_finite The colon is important If this command prints output similar to the following NOTICE halt_non_finite Integration halted then at least the check is not optimized away and halt_non finite is more likely to work 204 LANGUAGE REFERENCE When using GCC xmds enables the ffast math flag In GCC version 3 3 this seems to pose problems with halt non finite Removing the ffast math flag should solve the problem Alternatively keeping this flag and adding the fno unsafe math optimizations flag after the ffast math flag seems to work and should produce fas
134. f transverse dimensions here but in our example things are a bit simpler 1 2 A MORE COMPLEX SIMULATION 23 We now have to tell xmds the number of grid points of the lattice in this dimension and over what domain in this dimension the grid is defined To do these two things we use the lt lattice gt and lt domains gt tags In our simulation here we want to sample from x 1 to x 1 and use 100 points Therefore we set lattice to 100 and domains to 1 1 Notice that the lt domains gt tag is defined by using an ordered pair syntax xmds expects to see the domain for each transverse dimension defined as an ordered pair i e two comma separated values enclosed in parentheses if there is more than one transverse dimension then the domains for each are defined by a list of space separated ordered pairs The last thing we need to mention before discussing the lt vector gt tag is that the number of samples is set to the same value as that in Section 1 1 i e 1 and that the name of the field as per usual is main The lt field gt element now looks like lt field gt lt name gt main lt name gt lt dimensions gt x lt dimensions gt lt lattice gt 100 lt lattice gt lt domains gt 1 1 lt domains gt lt samples gt 1 lt samples gt lt more xmds code to cone s lt field gt We now need to specify the lt vector gt element This is much the same as previously discussed in Section but with some
135. faults to main Example simulation sequence integrate lt k operators vectors main vci lt vectors gt X k operators lt integrate gt lt sequence gt lt simulation gt 11 22 2 12 2 constant optional lt constant gt bool lt constant gt Contains boolean Subelement of Path to tag lt simulation gt lt sequence gt lt integrate gt lt k operators gt lt constant gt Description Tells xmds whether or not the k space operators are constant over the course of the simulation in other words they don t depend upon the propagation dimension If they are constant then giving a value of yes for this tag speeds up the simulation as more efficient code can be used Defaults to no Example 11 22 SEQUENCE 207 lt simulation gt lt sequence gt lt integrate gt lt k_operators gt lt constant gt yes lt constant gt lt k_operators gt lt integrate gt lt sequence gt lt simulation gt 11 22 2 12 3 operator names required operator names string variable Name string variableName operator names Contains array of strings Subelement of Path to tag simulation sequence integrate k operators operator names Description The names of the k space operators as they appear in the CDATA block within the k operators element This is a space separated list of strings of the operator names Example lt
136. fficient const double sigma 0 1 std dev of initial Gaussian const double x0 0 0 mean position of initial Gaussian ME lt globals gt Field to be integrated over gt ed lt name gt main lt name gt lt dimensions gt x lt dimensions gt lt lattice gt 100 lt lattice gt lt domains gt 1 1 lt domains gt lt samples gt 1 lt samples gt lt vector gt lt name gt main lt name gt lt type gt complex lt type gt lt components gt T lt components gt lt fourier_space gt no lt fourier_space gt lt CDATA T rcomplex exp x x0 x x0 2 0 signa sigma sigma sqrt 2 0 M_PI 20 0 JN lt vector gt 28 STARTING FROM SCRATCH lt field gt lt gt lhe siequence of integrations to perform gt lt sequence gt lt integrate gt lt algorithm gt RK4EX lt algorithm gt lt interval gt 1 lt interval gt lt lattice gt 1000 lt lattice gt lt samples gt 50 lt samples gt lt k_operators gt lt constant gt yes lt constant gt lt operator_names gt L lt operator_names gt lt CDATA L kappa kx kx gt lt k_operators gt lt CDATA dido LETE T8 lt integrate gt lt sequence gt lt U The OUEN to generate gt lt output gt lt filename gt diffusion xsil lt filename gt lt group gt lt sampling gt lt fourier_space gt no lt fourier_space gt lt lattice gt 50 lt lattice gt lt moments gt temperature lt moment
137. fomit frame pointer o advection advection cc I home cochrane bin lstdc 1m lxmds L home cochrane bin lfftw threads lfftw advection ready to execute then advection Performing fftw calculations Standing upon the shoulders of giants Importing wisdom Making forward plan Making backward plan Keeping accumulated wisdom Finished fftw calculations Beginning full step integration Sampled field for moment group 1 at t 0 000000e 00 Sampled field for moment group 1 at t ML 2 000000e 0 Sampled field for moment group 1 at t 4 000000e 02 lt snip gt Sampled field for moment group 1 at t 9 6O00000es gt 01 Sampled field for moment group 1 at t 9 800000e 01 Sampled field for moment group 1 at t 1 000000e 00 maximum step error in moment group 1 was 7 801465e 06 Time elapsed for simulation is O seconds That was fast eh And the error isn t too bad at about 10 so we can be vaguely confident of the results Lets look at them now in both Matlab or Octave and Scilab 3 4 1 Matlab and Octave Using xsil2graphics we generate the Matlab or Octave script by the command xsil2graphics advection xsil and then running the following commands in Matlab or Octave gt gt advection gt gt mesh t_1 x_1 amp_1 gt gt label 7 2 gt gt ylabel x gt gt zlabel A Doing all this should produce something very similar to that in Figure 62 USING A TEMPLATE
138. format binary assignment in the lt filename gt tag and then xmds would expect the data to be a string of double precison numbers in the same order as that given above 2 6 3 Importing complex data If we want to import complex data we just specify the real then imaginary parts sequentially as pairs of data Imagine that we now have two vectors so we don t have to consider so many vectors called x and y They have values of just for the sake of argument 1 2 2 0i 7 5 0 0i 5e 2 10i 7e10 8e 7i lt II and they will be organised in the input file as follows real x 0 imag x 0 real y 0 imag y 0 real x 1 imag x 1 real y 1 imag y 1 which is 1 2 2 0 5e 2 10 7 5 0 0 7e10 8e 7 With complex data the binary input method is slightly different The assignment to the format attribute is the same i e binary however instead of separating the real and imaginary parts of the complex numbers that are to be read in one just has the binary representation of the complex number to be read So in a sense the binary input of complex data is exactly the same as that of double data except that the data is complex and not double which seems obvious but it sort of had to be said 44 EXTRA AND ADVANCED FEATURES 2 7 Command line arguments Do you want to run your simulation many many times ranging over several different global parameters If the answer is yes then the command line argument f
139. g your code At the moment xmds version 1 5 1 there are six algorithms to choose from RK4EX RK4IP ARK45EX ARK451P SIEX SIIP RK4EX is a fourth order Runge Kutta in the explicit picture RK4IP is a fourth order Runge Kutta in the interaction picture the ARK45EX and ARK45IP are the corresponding adaptive time step Runge Kutta Fehlberg methods SIEX is the semi implicit method in the explicit picture and SIIP is the semi implicit method in the interaction picture For more information about the specifics of these algorithms and techniques see Section In solving our problem we ll use the fourth order Runge Kutta in the explicit picture because we aren t using any Fourier transforms and the explicit picture is fine for our purposes here Therefore we specify within the lt integrate gt tags the line lt algorithm gt RK4EX lt algorithm gt telling xmds what algorithm to use The next things xmds needs to know are the length of the integration interval the total number of steps to take and the number of samples for each output moment to take within these steps These items are denoted by the lt interval gt lt lattice gt and lt samples gt tags respectively The integration interval combined with the number of steps gives the step size internally used by xmds We ll choose some fairly arbitrary numbers here an interval length of 10 a large number of lattice points namely 10000 which should give us a nice small step size and we ll
140. gt lt CDATA const double kappa 0 1 diffusion coefficient const double sigma 0 1 std dev of initial Gaussian const double x0 0 0 mean pos of initial Gaussian 1 gt lt globals gt lt I remainder of dattiusion simulation xmds code gt lt simulation gt Notice that we ve commented out the kappa variable using the C line comment style This is just to remind us that kappa used to be there and is no longer and what it was when we originally wrote the simulation It can be a good idea to keep this kind of information around if you want but it isn t necessary and because it s a comment it will be ignored by the C C compiler Of course if you don t comment the global declaration out then the C C compiler will throw an error and your simulation won t compile Running xmds on the file diffusion xmds now gives a simulation binary that can accept arguments You can try it out by running the simulation like so diffusion kappa 0 2 where we have run the diffusion simulation with kappa now set to 0 2 xmds uses the GNU getopt set of functions to implement arguments and as such supports both short and long option names Therefore the above example could have been run as diffusion k 0 2 So at the simplest level xmds takes the long form of the argument name as the actual name of the variable and takes the first character of the variable name for the short form of the argument But what happens
141. gt lt RK4EX RK4IP SIEX SIIP gt lt iberations gt 3 lt itberations gt l default 3 for SI algs gt lt interval gt 10 lt interval gt S not Mens da main E gt lt lattice gt 1000 lt lattice gt no pointues alia mea Chan gt lt samples gt 100 lt samples gt lt no pts in output moment group gt lt CDATA dz_dt i omega z i sqrt 2 0 gam n 1 z e lt integrate gt lt sequence gt To generate the lt output gt block we choose the format to be ascii and remove the now superfluous precision attribute The sampling isn t done in Fourier space so we can set that to no we want to sample on a lattice of 100 points and take the moments realz and imagz which are the real and imaginary parts of z respectively These are then defined in 68 STOCHASTIC SIMULATIONS AND MPI the CDATA block by the following code lt CDATA realz real z imagz imag z 1 gt We therefore have an lt output gt block of lt I The output to generate gt lt output format ascii gt lt group gt lt sampling gt lt fourier_space gt no lt fourier_space gt lt sample in k space gt lt lattice gt 100 lt lattice gt AMS me pointes BOQ Sample c lt moments gt realz imagz lt moments gt lt names of moments gt lt CDATA realz real z imagz imag z 1 gt lt sampling gt lt group gt lt output gt and a final output simulation of lt xml versi
142. gt lt integrate gt no SIP integer defaults to 3 lt tolerance gt lt integrate gt yes ARK45 string lt max_iterations gt lt integrate gt no ARK45 integer defaults to infinity lt min_time_step gt lt integrate gt no ARRK45 89 double defaults to 1e 13 lt cutoff gt lt integrate gt no ARK45 string defaults to 1e 3 lt smallmemory gt lt integrate gt no ARK45IP yes no defaults to no lt halt_non_finite gt lt integrate gt no yes no defaults to no lt moments gt lt sampling gt yes array of strings post propagation default value arg yes depends upon type declaration benchmark simulation no yes no defaults to no use wisdom simulation no yes no defaults to no binary output simulation no yes no defaults to no deprecated use double simulation no yes no defaults to no deprecated use prefs gt simulation no yes no defaults to yes Table 8 3 The assignment elements 8 3 C coding within elements xmds is all about transplanting equations into the necessary multi dimensional loop struc ture to solve them which has the direct consequence that it is necessary for your equations to be written in standard C syntax To begin with variable names may not start with a num ber and may not be a reserved keyword such as void int long double complex for while if switch return etc This applies to the variable name
143. h that the additional memory use when used in a deterministic simulation with MPI is equal to the size of the field stored on any given node and only while the XSIL file is being written instead of the total field size for the entire simulation The syntax for a breakpoint element this should be in a sequence element is lt breakpoint gt filename XSIL filename for output e g simulation xsil gt lt filename gt lt fourier_space gt lt yes no gt lt fourier_space gt lt vectors US inst or vectors to be saved to the tale gt gt gt lt vectors gt lt breakpoint gt Using a template In this tutorial we re going to hack an xmds template to pieces to write a simulation script This is possibly one of the easiest ways to make an xmds simulation script almost everything has been done for you and all that is left for you to do is to translate the equations into something a C C compiler could understand This tutorial builds upon the skills learnt in Chapter 1 and from the knowledge of the existence of at least of the extra tags discussed in Chapter 3 1 The advection equation But first the problem Here we ll solve the one dimensional advection equation 0 3 400 05 A0 3 1 where A z t is the field to evolve according to the differential equation x is the spatial dimension t is time and v is the velocity of the wave As with the two examples discussed in Chapter 1
144. h zero The second step is to use define macros to map the equation variables to the internal data arrays within the generated code Thus for the example above using an explicit picture method the output code in the derivative calculation routine becomes 120 DEVELOPMENT AND PROGRAM STRUCTURE xmdsArg xmdsArgElement xmdsBreakPoint xmdsField xmdsMomentGroup xmdsintegrateARK45EX xmdsGlobals xmdsintegrateARK45 xmdsintegrateARK451P xmdsElement xmdsOutput l xmdsintegrateRK4 xmdsintegrateRK4EX xmdsSegment xmdsutility xmdsVector xmdsSimulation xmdsFilter xmdsintegrate K l xmdsintegrateEX xmdsintegrateRK4IP xmdsintegratelP xmdsintegrateSIEX xmdsintegrateSl xmdsintegrateSIIP xsilField xmdsVectorElement FIGURE 7 2 xmds class hierarchy void _sgi_calculate_delta_a const double amp _step complex dphi_dz _main_main_go_space 0 _main_sgi_coterms_go_space 0 unsigned long _main_main_pointer 0 unsigned long _main_sgi_coterms_pointer 0 double t main xminO for unsigned long _i0 0 iO main latticeO _i0 kxkkkkkkkkkk propagation code KKKKKKKKKKKK dphi_dz _sg1_coterm_L_phi i phi phi phi TR RRO ORO ok ook okokok ok ok ORR ORI OK R R K _main_main _main_main_pointer 0 dphi dz step main main pointer main main ncomponents main sgi coterms pointer main sg coterms ncomponents t main dxO0 with
145. have been put all through the document to give a better idea of what the tags do and what their default value is if any lt xml version 1 0 gt lt simulation gt I Global system parameters and functionality 2 name lt name gt lt the name of the simulation gt author lt author gt lt the author of the simulation gt description lt description Oo macs mie MesmmusiParmomWseSppoesedaoNdo gt lt description gt lt prop_dim gt lt prop_dim gt lt name of main propagation dim gt lt stochastic gt no lt stochastic gt lt defaults to no gt lt lt gt gt these three tags only necessary when stochastic as yes gt lt paths gt lt paths gt SL MO OF pathes gt lt seed gt 1 2 lt seed gt lt l seeds for rand mo gen gt lt noises gt lt noises gt ll mo of noises gt lt use_mpi gt no lt use_mpi gt x desamlts we mo gt lt error_check gt yes lt error_check gt lt l defaults wo yes gt lt use_wisdom gt yes lt use_wisdom gt lt der at sito mo gt lt benchmark gt yes lt benchmark gt lt I defaults to no gt lt use_prefs gt yes lt use_prefs gt e adeat MONSIG M 3 2 THE TEMPLATE CODE 53 EL glo lt C des lt gl lt lt fie lt n Global variables for the simulation gt bals gt DATA obals gt Eield to be integrated over gt gt ld ame main lt name gt dimensions
146. he length of the x domain SD ZA 5 4 x m uli 5 INN TONO Re o gt a p ps lt ER E 27 val N un Ji UNI 0 d bit l a t e A o 05 a Single path e 2 gt bee Hag P E 1 e 27 n e 10 12 2 Joel Z Pu m SAGs d Meee SE ER OR ee SRR ss eR BO V EMAN bs LIO Z D ie Wine ee 2 LY ae Z iy M7 71 7 A Hi 6 ao E b 1024 path mean FIGURE 10 3 Results for fibre xmds One important issue here is that the variance of the noise terms now scales with the product of the number of lattice points for a given domain Hence changing to a finer lattice actually increases the single trajectory error and the relationship between the error and the lattice product will depend on the order of the spatial derivatives The only way to overcome this is to reduce the integration step size which added to the fact that there are more lattice points in the first place means that fine lattice resolution in a stochastic PDE is computationally very expensive 10 4 tla xmds lt xml version 1 0 gt lt Two Level Atom Example simulation to crossi propagating et terdo m gt lt simulation gt lt prop_dim gt z lt prop_dim gt lt globals gt illustra
147. he simplified Kubo oscillator model taken from Gardiner 3 Z ilo VAE 4 1 This system is a model of an oscillator with a mean frequency w perturbed by a noise term t The parameter y describes the strength of the noise perturbation Examples of where this could be used to model an actual physical system are in the theories of magnetic resonance and laser physics and in single molecule spectroscopy 4 1 Without MPI First off let s solve this problem without MPI For those of you who don t know MPI is the Message Passing Interface and is the current de facto standard for performing computer simulations in parallel MPI is a very powerful library of routines that allow one to perform different parts of a simulation concurrently on different processors and be able to handle the relevant communication between processors or to run the same simulation with different parameters on different processors thereby in both situations speeding up the solution of the given problem Other similar systems exist such as PVM the Parallel Virtual Machine but MPI is effectively an extension of this system and is very well developed and works well on supercomputers and cluster systems MPI is also what xmds uses to parallelise its simulations so it gives us further motivation to only discuss it here We now take our template and hack around with it a bit Our propagation dimension 65 66 STOCHASTIC SIMULATIONS AND MPI is time so we set lt p
148. heck that for example the simulation isn t running off the grid or that the behaviour is wrong and the simulation needs to be terminated This way much time can be saved waiting for the result of a long simulation that needs to be re run anyway Another use of breakpoint elements is to save the state of some or all vectors to an XSIL file for loading by another simulation as discussed earlier in Section 2 6 1 of this chapter The naming convention for the vectors and components of complex vectors in the XSIL file produced is the same as that used for loading XSIL files as described in Section 2 6 1 Although creating an XSIL file to be used for initialising another simulation can be achieved almost as easily with an output moment group breakpoints should be used in stead of output moment groups for large deterministic simulations that use MPI Currently xmds 1 5 2 because of the way in which output moment groups are sampled with MPI each node allocates the entire memory required for sampling each output moment group This means that sampling the entire field for a large simulation that will not fit into the memory of a single node is impossible and hence creating an XSIL file from which the sim ulation could be continued is also impossible As the intended use of moment groups is that they should be used to sample a small amount of data an alternative solution was required for this situation Breakpoint elements have been designed suc
149. hen a template code will be written to either standard output i e the terminal or to file if a filename is given after the template flag For instance if one entered the following at the command line xmds template then you would see several lines of xmds code scroll past What s the use in it just scrolling past you say Well you could pipe this output to file like so 4 xmds template gt new xmds file xmds Nice eh However it is possible to save on keystrokes by getting xmds to make the file directly saving you from having to use shell related commands to save the output To save a template directly to file just enter this command xmds template new xmds file xmds and this will make a new file for you called new xmds file xmds in the same directory as xmds was called Why have the ability to just send the template to the screen Well doing so doesn t slow things down and it gives an extra level of flexibility and we re here to make your life as the user of xmds as easy as possible Stochastic simulations and MPI One of the most powerful features of xmds is its ability to automatically parallelise stochastic funky way to say random simulations To illustrate this power we re going to have a look at an example of a stochastic differential equation SDE both with and without the use of MPI the Message Passing Interface used for running parallel simulations The physical model we ll be investigating is t
150. her the Program or any derivative work under copyright law that is to say a work containing the Program or a portion of it either verbatim or with modifications and or translated into another language Hereinafter translation is included without limitation in the term modification Each licensee is addressed as you Activities other than copying distribution and modification are not covered by this Li cense they are outside its scope The act of running the Program is not restricted and the output from the Program is covered only if its contents constitute a work based on the Program independent of having been made by running the Program Whether that is true depends on what the Program does 1 You may copy and distribute verbatim copies of the Program s source code as you receive it in any medium provided that you conspicuously and appropriately publish on each copy an appropriate copyright notice and disclaimer of warranty keep intact all the notices that refer to this License and to the absence of any warranty and give any other recipients of the Program a copy of this License along with the Program You may charge a fee for the physical act of transferring a copy and you may at your option offer warranty protection in exchange for a fee 2 You may modify your copy or copies of the Program or any portion of it thus forming a work based on the Program and copy and distribute such modifications or work under the terms of Sectio
151. here However preferences are on by default anyway and switching them off probably won t be of much use so we ll get rid of the lt use_prefs gt tag which will set it to the default value of yes We ve seen in Chapter 2 how one can use command line arguments with the lt argv gt tagset At this stage of simulation writing it s a good idea just to get things going and not worry about varying these parameters just yet so what we ll do is set them in the lt globals gt block and add the lt argv gt stuff in later if we want to The first chunk of code is now I Global system parameters and functionality gt lt name gt advection lt name gt lt the name of the simulation gt lt author gt Paul Cochrane lt author gt lt the author of sim gt lt description gt Solves the one dimensional advection equation for an initial cosine modulated Gaussian pulse Adapted from A L Garcia Numerical Methods in Physics 1994 lt description gt lt prop_dim gt t lt prop_dim gt lt name of main propagation dim gt lt error_check gt yes lt error_check gt s defaults Lo e ME use wisdom yes use wisdom lt U dsremiles 0 mo M s lt benchmark gt yes lt benchmark gt Se AGS wo mO 3 3 2 Global variables for the simulation Our global variables are v k zo and o from inspection of the equations in Section 3 1 therefore we need to use the following code for the lt globals gt block
152. his implies that the lt paths gt lt seed gt and lt noises gt tags are unnecessary Also we don t need the lt use_mpi gt tag either because we ve not met it before the simulation is not stochastic so making it use MPI the message passing interface for parallel simulations is a waste of time and the default setting is off so again we rip this tag out 3 3 RIPPING IT TO BITS 55 We might as well leave error checking on it s the default anyway and leave the lt error_check gt tag the way it is This is especially important at the early stages of writing a simulation as one can then get an indication of the discretisation error of the simulation and this in formation can tell us how to tweak the simulation parameters if necessary Leaving the lt error_check gt tag in the code is a good idea as we ll probably want to switch it off at some later stage and if we leave it in the code we can remember to do this if we want We definitely want to use FFTW s wisdom feature This is because we have a differential in space and it is nice and easy to define the operator in k space implying that we ll need to use Fourier transforms and using wisdom speeds up the startup time of our simulations Hence we leave lt use_wisdom gt as is Benchmarking the code doesn t take up much room in the script nor does it take up much time in the code and it can be interesting to have around so we ll leave the benchmark tag t
153. i 4 i 0 5 phi 2 phi 2 ixgam13 gamloss132 densixphi 4 dphi_dt 5 L2p phi 5 i Vr 4 gamlossi1 phi 5 i gameff gamloss12 gdensi ichippy ixgami3 gamloss132 gdens3x phi 5 i conj phi 5 phi 6 dphi_dt 6 L4p phi 6 i Vr 5 gamloss31 phi 6 i gam33 gamloss32 gdensS3 phi 6 i 0 5 phi 5 phi 5 i gami3 gamlossi32 gdensi phi 6 gt lt integrate gt lt sequence gt lt he output to generate gt gt gt 172 MORE EXAMPLES lt output format binary precision double gt lt group gt lt sampling gt Tour ie nts pace no no no lt fourier_space gt lt lattice gt 16 1 1 lt lattice gt lt moments gt atoms molecules gatoms gmolecules lt moments gt lt CDATA atoms phi 2 phi 1 molecules phi 4 phi 3 gatoms conj phi 5 phi 5 gmolecules conj phi 6 phi 6 I gt lt sampling gt lt group gt lt group gt lt sampling gt lt fourier_space gt no no no lt fourier_space gt lt lattice gt 0 16 0 lt lattice gt lt moments gt rn_1 rn_2 grn_1 grn_2 excitedn lt moments gt lt CDATA rn_1 phi 2 pa1 1 rn_2 phi 4 phi 3 grn_1 conj phi 5 phi 5 grn_2 conj phi 6 phi 6 excitedn g g 4 phi 2 phi 2 phi 1 phi 1 F F phi 4 phi 3 Fxg 2 phi 2 phi 2 phi 3 phi 1 phi 1 phi 4 1 gt lt sampling gt lt group gt lt group gt lt sampling gt lt fourier_space gt no no no lt fourier_space g
154. i3 Uoh13 chi const double gam33 Uoh33 chi const double gameff Uoh11 biggamma chi const double gamlossii loss11 2 chi const double gamloss12 loss12 chi const double gamloss31 10ss31 2 chi const double gamloss32 10ss32 chi const double gamloss132 loss132 chi const double cnoise noise sqrt 2 0 1 gt lt globals gt lt argv gt lt arg gt lt name gt kjoek lt name gt lt type gt double lt type gt lt default_value gt 1 0e6 lt default_value gt 10 6 HIGHDIM_VECTOR_VERSION XMDS 169 lt arg gt lt arg gt lt name gt joekappamax lt name gt lt type gt double lt type gt lt default_value gt 1 0e2 lt default_value gt lt arg gt lt argv gt Field to be integrated over gt lt field gt lt dimensions gt x y z lt dimensions gt lt lattice gt 16 16 16 lt lattice gt lt domains gt 1 2e 4 1 2e 4 1 2e 4 1 2e 4 8 0e 3 8 0e 3 lt domains gt lt samples gt 1 1 1 lt samples gt lt vector gt lt name gt vcl lt name gt lt type gt double lt type gt lt components gt Vr 5 lt components gt lt fourier_space gt no no no lt fourier_space gt lt CDATA Vr 1 omegax romegax x x omegay omegay y y omegazromegaz zx z Vr 2 0 5 M Vr 1 hbar chi gameff gami3 2 2 dx dy dz Vr 3 M Vr 1 hbar chi gami3 2 gam33 2 dx dy dz Vr 4 0 5 M Vr 1 hbar chi Vr 5 M Vr 1 hbar chi gt lt vector gt lt vector gt lt name gt main lt n
155. iables will then be mentioned in the CDATA block following this tag Example lt simulation gt lt output gt lt group gt lt sampling gt lt moments gt pow_dens lt moments gt lt sampling gt lt group gt lt output gt lt simulation gt 11 23 2 2 post_propagation optional lt post_propagation gt xmds tags lt post_propagation gt Contains ass CDATA Subelement of Path to tag lt simulation gt lt output gt lt group gt lt post_propagation gt Description Container for tags describing any post propagation processing of the data that should be done prior to output to file Example lt simulation gt output group post propagation S xmds tags gt lt post_propagation gt lt group gt lt output gt lt simulation gt 11 23 OUTPUT 221 11 23 2 2 1 fourier space post propagation required lt fourier_space gt bool bool lt fourier space Contains array of booleans Subelement of post propagation Path to tag simulation output group post propagation lt fourier_space gt Description Whether or not the post propagation is performed in Fourier space This is a list of yes no entries for the propagation dimension and as many remaining transverse dimensions Example simulation output group post propagation X fourier space no fourier space X post propagation lt group gt
156. id and will automatically terminate your rights under this License How ever parties who have received copies or rights from you under this License will not have their licenses terminated so long as such parties remain in full compliance 5 You are not required to accept this License since you have not signed it However nothing else grants you permission to modify or distribute the Program or its derivative works These actions are prohibited by law if you do not accept this License Therefore by modifying or distributing the Program or any work based on the Program you indicate your acceptance of this License to do so and all its terms and conditions for copying distributing or modifying the Program or works based on it 6 Each time you redistribute the Program or any work based on the Program the recipient automatically receives a license from the original licensor to copy distribute or modify the Program subject to these terms and conditions You may not impose any further restrictions on the recipients exercise of the rights granted herein You are not responsible for enforcing compliance by third parties to this License 7 If as a consequence of a court judgment or allegation of patent infringement or for any other reason not limited to patent issues conditions are imposed on you whether by court order agreement or otherwise that contradict the conditions of this License they do not excuse you from the conditions of thi
157. in the home stretch now All we need to do is tell xmds what to spit out at the end of the simulation We ll use ascii output a bit more portable and remove the precision assign ment we won t sample the output moments in Fourier space we only have one dimension to sample namely x and we ll use 50 points here there is only one moment to output and that is the field amplitude which we ll call amp and finally the code to sample the output moment is simply equal to the amplitude of the field which is just A so the code put into the CDATA block is amp A The lt output gt element code is therefore lt Ui The output to generate gt lt output format binary precision single gt lt group gt lt sampling gt lt fourier_space gt no lt fourier_space gt lt sample in k space gt lt lattice gt 50 lt lattice gt SL 10 Pus BO sample gt lt moments gt amp lt moments gt lt names of moments i lt CDATA amp A qus lt sampling gt lt group gt lt output gt 3 3 RIPPING IT TO BITS 59 3 3 6 The final program And that s it The finished simulation is this Global system parameters and functionality gt lt name gt advection lt name gt lt the name of the simulation gt lt author gt Paul Cochrane lt author gt lt the author of sim gt lt description gt Solves the one dimensional advection equation for an initial cosine modulated Gaussian pulse Adapted f
158. in the very general form 9 a D s Oa Od These come under three sub classifications 5 according the the parameters A B and C 1 The PDE is said to be elliptic in regions of the z y plane where AC B 0 2 The PDE is said to be parabolic in regions of the z y plane where AC B 0 3 The PDE is said to be hyperbolic in regions of the z y plane where AC B 0 This model for classifying PDEs is only two dimensional but many problems in physics exhibit up to four dimensions However the majority of PDEs in physics are able to be expressed in a form where one of the dimensions is special special in that there are only ever first order derivatives involving that dimension A two dimensional form of such a PDE looks like Con sat s f ena tes 62 and such an equation falls into the category of parabolic by the above classification For this reason many three and four dimensional PDEs in physics are referred to as parabolic 81 82 NUMERICAL MODELLING THEORY although technically this classification exists only for PDEs in two dimensions The more general form for such parabolic like differential equations is o Ox i x f Pe a x 6 3 where the functionals ft may include first or second order partial derivatives of the transverse dimensions z 7 This set of equations can be generalised to allow the functionals f to include partial derivatives of any order Given that th
159. ing n11 bec11 recon testing n12 bec12 recon testing n13 bec13 recon testing n14 bec14 recon testing n15 bec15 Woo hoo recon has completed successfully This means that you will most likely be able to boot LAM successfully with the lamboot command but this is not a guarantee See the lamboot 1 manual page for more information on the lamboot command If you have problems booting LAM with lamboot even though recon worked successfully enable the d option to lamboot to examine each step of lamboot and see what fails Most situations where recon succeeds and lamboot fails have to do with the hboot 1 command that lamboot invokes on each host in the hostfile If everything has run successfully then you can start MPI on each of the nodes other computers by using the lamboot command This gets each of the other computers ready to get input from a parallel job and doesn t actually start doing any computation You should see output similar to this lamboot v lamhosts LAM 6 5 6 MPI 2 C ROMIO University of Notre Dame Executing hboot on 10 bec00 1 CPU Executing hboot on ni becO1 1 CPU Executing hboot on n2 bec02 1 CPU Executing hboot on n3 bec03 1 CPU Executing hboot on n4 bec04 1 CPU Executing hboot on n5 bec05 1 CPU Executing hboot on n6 bec06 1 CPU Executing hboot on n7 becO7 1 CPU Executing hboot on n8 bec08 1 CPU Executing h
160. ing one to do some complex tasks on the data sampled in the running of the simulation this is put after the lt sampling gt tag which is coming up however this is more involved than the discussion here so we won t be using it A good place to look if you re at all interested is the examples directory in the distribution or see the script repository on the smds web page After looking for the lt group gt tag xmds then expects to see a lt sampling gt tag within that which defines how the group is to be sampled This is also just a container for more specific tags Just as an aside although it may seem a pain at this stage to have so many containers for other containers and tags and so on this gives the entire document a nice structure where one builds up a simulation from nice bite sized pun not intended chunks Also one really doesn t go to the pain of writing a simulation from scratch very often so this shouldn t be a big issue when you finally get to writing other and more complex and more interesting simulations Within the lt sampling gt tag xmds expects to see a lt fourier_space gt tag for each transverse dimension in the simulation which tells xmds whether or not to Fourier transform the dimension before being sampled and written out to disk For our example here we don t have any transverse dimensions so we won t use this tag Now we need to tell xmds the names of output moments we want to calculate and the C code it
161. inning full step paths Starting path 1 Starting path 2 lt snip gt Starting path 1023999 Starting path 1024000 maximum step error in moment group 1 means was 9 161489e 05 Time elapsed for simulation is 1479 seconds Hmmm that took a bit longer didn t it That s just over 24 minutes Surely we can do better than that This is where we need to use many computers to solve the problem and so this is where MPI comes in 4 2 With MPI Now we hope that by farming each path off to a different CPU that we can get an improve ment in speed Let s try it and find out At this stage we should expect that since each individual path doesn t take very long to run then the main overhead is going to be communication between farming a job off and sending the results back In a more complex simulation it would be normal and in fact 74 STOCHASTIC SIMULATIONS AND MPI much better for each node to spend quite some time working on an individual path before sending the information back such a scenario is far more efficient computational use of CPU time In any case we should see a shorter wall time the amount of time taken as far as the clock on the wall is concerned for the job to complete To implement the use of MPI make sure that your copy of xmds is built with MPI enabled To do this one merely needs to run the configure script with the enable mpi option specified From that point the only thing you need to change within your XMDS s
162. ion gt or conversely just before the lt globals gt section in the situation considered here these locations are one in the same however this is not the case in general The script is now lt xml version 1 0 gt lt simulation gt lt name gt lorenz lt name gt lt author gt Paul Cochrane lt author gt lt description gt Lorenz attractor example simulation Adapted from the example in Numerical methods for physics by Alejandro L Garcia page 78 ist ed lt description gt lt prop_dim gt t lt prop_dim gt lt gt yet more xmdis code eto come gt gt gt lt simulation gt This problem we are trying to solve has several constants namely o r b zo yo and zo These variables are going to be used again and again in the simulation therefore it makes sense to put them into the next element necessary to describe a simulation in xmds the lt globals gt tag We specify these constants using C C syntax in an XML CDATA block which is enclosed within the lt globals gt element The xmds code for this is lt globals gt lt CDATA const double sigma 10 0 const double b 8 0 3 0 const double r 28 0 const double xo 1 0 initial conditions const double yo 1 0 initial conditions const double zo 20 0 initial conditions gt lt globals gt 1 1 A BASIC SIMULATION 7 There are a couple of points to note here Firstly and most importantly the syntax of the code within
163. ired accuracy many times faster than the semi implicit method 6 4 5 The Fourth Fifth Order adaptive Runge Kutta Method in the Explicit Picture Once more we use the Fourier Transform or Spectral method to determine the transverse derivatives In order to minimize memory usage the current contents of the arrays containing the final solution are reused to calculate some of the k vectors The optimized procedure after Section 6 3 5 then is 1 Calculate if required 2 ax a 3 ay a 4 p F L 1 k F ax 5 ax AN 2 ak p 6 a a ax 7 a a cax 8 x x ash 92 NUMERICAL MODELLING THEORY 9 ax a bajar 10 p F L 2 k1 F ax 11 ax AN 2 ax p 12 a 1 b31 c1 a1 b31 c1a 032a 13 z x a3 az h Note co ch 0 14 p F EOS k i Z akx 15 a A N zapis 16 az 1 d41 c1 ar b41 c1a 042a bazas 17 a a cza 18 a a Gay 19 z x a4 a3 h 20 p FH D e s az 21 a hN x ar p 22 a a Car 23 a a car 24 ar 1 b5 c1 ar bsi cia bs2aK bsa bs103 c1 as bsa b51 4 c1 az 25 x x as a4 h 26 p Et DE Gr ki e az 27 a hN x ar p 28 a a car Note cs 0 29 ar fiar fax faa faa fsar fea 30 x zx ag as h 31 p PAUL qe mel 32 az hN x az p 6 4 NUMERICAL METHODS FOR PARTIAL DIFFERENTIAL EQUATIONS 93 33 a
164. is aya h 75 a h FIN 2 F7 aj 76 x x aig a15 h TT aj F IN x F a T8 ara a 79 a a cia cgay coa Cioad 0119 C128 f C13ag C14Ap C158j C168 80 aa Ame C aa c a Cac C198q Cie 015385 Cigay Cigan clsa Ci a The method requires 11 copies of the field and has a variable time step which is stochastically safe like its explicit picture counterpart Once again no rotations are required for steps 15 and 16 since they both occur at the end of the time step 6 4 14 Adding a Cross Vector The nonlinear functionals of Equations and may also include variables that are a function of the same transverse space but are which are governed by ODEs in one of the transverse dimensions These components need never be transformed to Fourier space and so ought to be separated from the main vector components for efficiency reasons In Equation they are referred to as the vector b which propagates in dimension z This vector has a definite value once the main vector a has a definite value This means that for most integration algorithms a current value for the cross vector may be calculated immediately prior to calculating the M functionals In the semi implicit method in the explicit picture the cross vector b is calculated using the semi implicit method and using linear interpolation in the x dimension to create the mid lattice point values for a that are required Simil
165. is used initially as the mid point value am a tm tm t h 2 and then the step 1 am a h flim aml 6 18 is iterated until convergence is obtained The next value for a is then An 1 2am s 6 19 Similar to the improved Euler method the semi implicit method is also a second order method It has the added advantage that a na ve inclusion of the noise terms in the function f produces an algorithm that converges to the Stratonovich integral with global error of the order of h 6 3 4 The Fourth Order Runge Kutta Method The fourth order Runge Kutta method has been the workhorse for numerical modelling for many years The method uses four evaluations of the derivative to account for the first second third and fourth derivatives making it a fourth order method The method is simply ki A f tn anl 6 20 k hf Sn bus a e ski 6 21 ks hf Sta et E ak ska 6 22 ka A f t i an ks 6 23 ME T ku ola Joke 6 24 6 25 Higher order methods are possible but usually the extra accuracy does not warrant the extra computation In other words this method usually achieves the fastest evaluation for a given propagation accuracy Note that even an optimised implementation of this algorithm will require about three times as much memory as a lower order method 86 NUMERICAL MODELLING THEORY Although empirically there is a subclass of Stratonovich stochastic problems that converge with
166. kpoint gt lt filename gt Description Sets the XSIL file to which the breakpoint information should be saved Example lt simulation gt lt sequence gt lt breakpoint gt lt filename gt blah xsil lt filename gt lt breakpoint gt lt sequence gt lt simulation gt 11 22 4 2 fourier_space breakpoint optional lt fourier_space gt boolean boolean lt fourier_space gt Contains array of booleans Subelement of Path to tag lt simulation gt lt sequence gt lt breakpoint gt lt fourier_space gt Description Tells xmds in which space the breakpoint is to be written in This is a list of yes no options for each vector Example 11 23 OUTPUT 215 lt simulation gt lt sequence gt lt breakpoint gt lt fourier_space gt yes no lt fourier_space gt lt breakpoint gt lt sequence gt lt simulation gt 11 22 4 3 vectors breakpoint required lt vectors gt string vectorName string vectorName lt vectors gt Contains array of strings Subelement of Path to tag lt simulation gt lt sequence gt lt breakpoint gt lt vectors gt Description The names of the vectors that will be saved to the XSIL file Example lt simulation gt lt sequence gt lt breakpoint gt lt vectors gt main vcl lt vectors gt lt breakpoint gt lt sequence gt lt simulation gt 11 23 output required lt output gt xmds tags lt output gt Contains
167. lation yet the list of instructions necessary to get a computer to calculate the solution is vast in comparison Going from the former to the latter is non trivial It can be done in a single step as far as the user is concerned with high level programming languages but most computer programmers know from experience that the numerical solution can be calculated more efficiently by writing comparatively more code in a low level script than less code in a high level script This is because the programmer knows the purpose of the program and can optimise it accordingly as the list of instructions is expanded On the other hand a compiler can only perform optimisations that it knows are guaranteed to work regardless of the purpose of the program By restricting the scope of programs to those with a known range of purposes one is able to write a compiler that generates more efficient output code In other words the possible input code is constrained to a very small set of instructions which even though they may be complex have a well defined expansion into a simpler instruction set Further it is not necessary to define the expansion down to the level of assembly code one only needs to define the expansion to a point where it no longer relies critically on the compiler s optimisation ability This is how xmds works transform from specialised input script to efficient C code and then let an ordinary compiler do the rest In choosing the high le
168. lation name which we ll call diffusion and we ll put the author name and a brief description of what the simulation is supposed to do The global variables we have in the problem we are solving are the diffusion coefficient k the standard deviation of the initial Gaussian distribution o and the positon of the mean of the Gaussian distribution zo These variables we ll call respectively kappa sigma and x0 The code for this looks like lt xml version 1 0 gt lt simulation gt lt Global system parameters and functionality gt lt name gt diffusion lt name gt lt author gt Paul Cochrane lt author gt lt description gt Solves the one dimensional diffusion equation for an initial Gaussian pulse Adapted from A L Garcia Numerical Methods in Physics 1994 lt description gt lt gt Global viarntables tor the simulation gt gt gt lt globals gt lt CDATA const double kappa 0 1 const double sigma 0 1 const double x0 0 0 113 lt globals gt lt more xmds code to come gt lt simulation gt which we put into a file called diffusion xmds 1 2 3 Describing the field Again the next thing to tell xmds about is the lt field gt element This time we have a transverse dimension which is in the x direction and we mention this using the lt dimensions gt tag like so lt dimensions gt x lt dimensions gt In the general case xmds expects a space separated list o
169. lement of Path to tag lt simulation gt lt sequence gt lt integrate gt lt samples gt Description The number of samples to take of each output moment group This is a space separated list of integers the number of which must be equal to the number of output moment groups Each integer must be a factor of the lt lattice gt assignment or 0 If set to O then the given moment group is not sampled Example lt simulation gt lt sequence gt lt integrate gt lt samples gt 50 lt samples gt lt integrate gt lt sequence gt lt simulation gt 11 22 2 12 k_operators optional lt k_operators gt xmds tags lt k_operators gt Contains operator nanes gt CDATA Subelement of Path to tag lt simulation gt lt sequence gt lt integrate gt lt k_operators gt Description Container for tags describing the k space i e Fourier space operators Example lt simulation gt lt sequence gt lt integrate gt lt k_operators gt lt U ssc vase gt lt k_operators gt 206 LANGUAGE REFERENCE lt integrate gt lt sequence gt lt simulation gt 11 22 2 12 1 vectors k_operators optional lt vectors gt string variableName string variableName lt vectors gt Contains array of strings Subelement of Path to tag simulation sequence integrate lt k operators vectors Description Vectors to be referred to in CDATA block De
170. ll refer to this script as lorenz xmds Next the xmds script needs a name Well to be honest it doesn t really need a name but it s a good idea to add the lt name gt tag if only for completeness If the lt name gt tag is not included xmds uses the name of the xmds script file but with the xmds extension removed as the name of the simulation so in this simulation if we were to not specify the lt name gt then the simulation name used by xmds would be lorenz However it is still a good idea to specify the name inside the simulation text and we do this here setting lt name gt to lorenz The script now becomes lt xml version 1 0 gt lt simulation gt lt name gt lorenz lt name gt lt gt yet more xmds code to come gt gt gt lt s m la b lons Notice that we have indented the lt name gt tag from the rest of the tags It is a good idea to indent tags that are nested within other tags so that one gets an idea of the document structure just from looking at the source code The next two tags that I recommend you add to your scripts are the lt author gt and lt description gt tags These tags aren t necessary for running a simulation or for actually calculating anything however they are good for documenting the simulation and providing extra information that may be helpful to others if you show your scripts to other people Also this information may actually be used in future versions of xmds to provide extra fun
171. lly better use a low order polynomial approximation 9 7 The filter element The lt filter gt element though not covered in this example is relatively straightforward Here is an example which retains a Gaussian profile of the field centered on t 0 lt filter gt lt vectors gt main lt vectors gt lt fourier_space gt no lt fourier_space gt lt CDATA phi exp t t 9 8 THE OUTPUT ELEMENT 145 11 gt lt filter gt As in the lt vector gt element the space in which the filter is to be applied is specified in a lt fourier_space gt assignment lt moment_group gt and lt functions gt elements may be used in the lt filter gt element just as they are used in the lt integrate gt element 9 8 The output element lt output gt lt filename gt nlse xsil lt filename gt lt group gt lt sampling gt lt fourier_space gt no lt fourier_space gt lt lattice gt 50 lt lattice gt lt moments gt pow_dens lt moments gt lt CDATA pow_dens phi phi BE lt sampling gt lt group gt lt output gt The lt output gt element is compulsory and defines what will be used as the output The lt filename gt assignment is optional From version 1 2 xmds can produce binary output as well as the default ascii format Both output formats will produce an ascii file in XSIL format 2 The default name for this file is the lt simulation gt lt name gt appended with xsil With ascii outpu
172. lt group gt lt sampling gt lt type gt complex lt type gt lt sampling gt lt group gt lt output gt lt simulation gt 11 23 2 1 2 fourier_space sampling required lt fourier_space gt bool bool lt fourier_space gt Contains array of booleans Subelement of Path to tag lt simulation gt lt output gt lt group gt lt sampling gt lt fourier_space gt Description A boolean telling xmds whether or not to sample in Fourier space This should be a space separated list of booleans one for each transverse dimension Example lt simulation gt lt output gt lt group gt lt sampling gt lt fourier_space gt no lt fourier_space gt lt sampling gt lt group gt lt output gt lt simulation gt 11 23 2 1 3 vectors sampling optional lt vectors gt string variableName string variableName lt vectors gt Contains array of strings Subelement of Path to tag lt simulation gt lt output gt lt group gt lt sampling gt lt vectors gt 11 23 OUTPUT 219 Description Space separated list of strings giving the names of the vectors to sample Defaults to main Example lt simulation gt lt output gt lt group gt lt sampling gt lt vectors gt main vcl lt vectors gt lt sampling gt lt group gt lt output gt lt simulation gt 11 23 2 1 4 lattice sampling optional lt lattice gt int int lt lattice gt
173. lt output gt group Description Container for tags describing the relevant moment group There must be at least one group element given Example 11 23 OUTPUT 217 lt simulation gt lt output gt lt group gt S mee TERE 2 lt group gt lt output gt lt simulation gt 11 23 2 1 sampling required lt sampling gt xmds tags lt sampling gt Contains moments CDATA Subelement of Path to tag lt simulation gt lt output gt lt group gt lt sampling gt Description Container for tags describing how to sample the moment group Example lt simulation gt lt output gt lt group gt lt sampling gt lt SU gt ames tags gt lt sampling gt lt group gt lt output gt lt simulation gt 11 23 2 1 1 type sampling optional lt type gt string of complex or double lt type gt Contains string interpreted as either complex or double type Subelement of Path to tag simulation output group sampling type Description The data type of the output data Defaults to complex if the vector that the output data depends on is of type complex and double otherwise Note that a type of double cannot be used with a post propagation tag as the fourier trans forms available to the post propagation tag require the output data to be of type complex 218 LANGUAGE REFERENCE Example lt simulation gt lt output gt
174. lty was developing 119 M Process command line arguments M Parse input file Write includes Y z Y Process simulation element Write defines Y Generate output code lt Write globals M Y x Call system compiler Write prototypes M Write routines FIGURE 7 1 xmds main procedure a C code style syntax which would allow a large range of complex PDEs to be encoded and in a manner that remains independent of the exact integrate algorithm chosen This was implemented using a two step process To begin with the user is required to write equations with any differential operators replaced by an operator component representation For example the NLSE shown in Equation 7 1 do Pe 9 would be written in C code style as dphi_dz L phi i phi phi phi where the operator L is defined in Fourier space as being L i kt kt 2 The first step is to search the user s equations for such operator component expressions and replace them with something else depending on the algorithm chosen In the explicit picture they are replaced with a reference to an internally calculated vector which is exactly the derivative calculated with the Fourier transform method whereas in the interaction picture the field is evolved in Fourier space according to the operator component expressions found and in the original code the text strings for these expressions are replaced wit
175. ly transplants these equations into the output C code which it then compiles and yet it had to be capable of representing equations of the form shown in Equation 8 3 Further it was desired that the C code form of such an equation must be the same regardless of the type of integration algorithm chosen This difficulty was compounded by the desire to include split operator algorithms which process the linear and nonlinear terms entirely separately The result that we have achieved here is quite elegant though there are a few caveats with regard to writing the equations The second area of difficulty was enabling the user to define exactly what they wanted as output The simplest solution would have been to save the raw field data to file as the simulation progressed and the user left to process it afterward in their own preferred plotting and calculation package However with stochastic problems the user usually wishes to calculate the average of some property or moment of the field over many different trajectories or integration runs Therefore it was necessary to enable to user to define such moments so that they may be calculated and averaged before the output data is written to file The source code for xmds is written in object oriented C otherwise known as C The main routine is very simple and its procedure is as shown on the left side in Figure The right side is performed by the simulation data class A non validating XML parser proce
176. m RK4EX lt algorithm gt lt RK4EX RK4IP ARK45EX ARKA45IP SIEX SIIP lt interval gt 1 lt interval gt O rar An maim Gaine gt lt lattice gt 500 lt lattice gt lt l no poinre ain mean Claw gt samples 50 lt samples gt lt no pts in output moment group gt lt k_operators gt lt constant gt yes lt constant gt yes no gt lt operator_names gt L lt operator_names gt lt CDATA L rcomplex 0 0 v kx 1 gt lt k_operators gt lt CDATA dA dt L A gt lt integrate gt lt sequence gt lt i The output to generate gt output format ascii gt lt group gt lt sampling gt lt fourier_space gt no lt fourier_space gt lt sample in k space gt lt lattice gt 50 lt lattice gt AO OBS WO Sewjyll a S lt moments gt amp lt moments gt lt gt names of moments gt gt gt lt CDATA amp A 11 lt sampling gt lt group gt lt output gt lt simulation gt Here is a link to the finished gzipped script file advection xmds gz on the xmds web site http www xmds org 3 4 MAKING THE SIMULATION AND GETTING RESULTS 61 3 4 Making the simulation and getting results As per usual we just need to run xmds on the simulation script and then run the simulation So here we go xmds advection xmds Output file name defaulting to advection xsil compiling g pthread 03 ffast math funroll all loops
177. m to Fourier space in the propagation dimension normal space in the x dimension and Fourier space in the z dimension when evaluating the post processing moments Finally there is another lt moments gt assignment with new names and more code in the content for processing the moments already derived at the sampling level As with most other code blocks functions of dimensions may also be included but only of the dimensions that haven t been collapsed plus the propagation dimension Finally it is only the real part of the lt post_propagation gt moments that are written as output or added to the stochastic sums before they are processed 9 8 THE OUTPUT ELEMENT 147 In the absence of a lt post_propagation gt element the real parts of the lt sampling gt mo ments are used instead 148 WORKED EXAMPLE NLSE XMDS 10 More Examples 10 1 ndparamp xmds S xml versions 11 0 v A NO De nene asma amet Amplifier lt Simulton formation for logical switching gt lt simulation gt lt name gt ndparamp lt name gt lt prop_dim gt z lt prop_dim gt lt error_check gt yes lt error_check gt lt globals gt lt CDATA const double el 350 const double e2 350 const double ri 1 const double r2 1 const double vyl 0 5 const double vy2 0 5 const double yci 0 2 const double yc2 0 2 const double tcl 0 const double tc2 0 double amp1 sqrt e1 2 M_P1 r1 r1 double amp2
178. m0 Pole o semple 2 lt moments gt realz imagz lt moments gt lt names of moments gt lt CDATA realz real z imagz imag z 11 lt sampling gt lt group gt lt output gt lt simulation gt The gzipped script code can be downloaded from the xmds web site http www xmds org via the link kubo tutorial xmds gz 4 1 1 Making the simulation and getting results You should be pretty good at doing this now so we ll just show you the output of the various operations and not explain much xmds kubo_tutorial xmds Output file name defaulting to kubo_tutorial xsil compiling g pthread 03 ffast math funroll all loops fomit frame pointer o kubo_tutorial kubo_tutorial cc I home cochrane bin lstdc 1m lxmds L home cochrane bin lfftw_threads lfftw kubo_tutorial ready to execute then kubo tutorial Beginning full step paths Starting path 1 Starting path lt snip gt Starting path 1023 Starting path 1024 maximum step error in moment group 1 means was 1 010483e 04 Time elapsed for simulation is 5 seconds That was pretty good for doing 1024 paths twice at different time steps eh Admittedly we re not pushing things much here We ll do so soon though 4 1 WITHOUT MPI 71 4 1 1 1 Matlab and Octave Using xsil2graphics we generate the Matlab or Octave script by the command xsil2graphics kubo tutorial xsil Running the following sequence of commands in Matlab or Octave
179. many maths libraries hence we use gam here double gam 0 1 perturbation strength double zo 1 initial value of z gt lt globals gt Field to be integrated over gt lt field gt lt name gt main lt name gt lt samples gt 1 lt samples gt ASES np e EST p OO Ene lt vector gt lt name gt main lt name gt lt type gt complex lt type gt lt I clima type ofi vectori gt lt components gt z lt components gt lt names of components gt lt fourier_space gt no lt fourier_space gt lt defined in k space gt lt CDATA z ZO Js lt vector gt lt field gt Su ihe seguencel of integrations to perform gt lt sequence gt lt integrate gt lt algorithm gt SIEX lt algorithm gt lt RK4EX RK4IP SIEX SIIP gt lt iterations gt 3 lt iterations gt lt default 3 for SI algs gt lt interval gt 10 lt interval gt AS ONE eve a mann came gt lt lattice gt 1000 lt lattice gt A AC poire sid mem Chim SSS lt samples gt 100 lt samples gt lt no pts in output moment group gt lt CDATA dz dt i omega z i sqrt 2 0 gam n 1 z 11 lt integrate gt lt sequence gt 70 STOCHASTIC SIMULATIONS AND MPI The out pu sto senerate gt gt gt lt output format ascii gt lt group gt lt sampling gt lt fourier_space gt no lt fourier_space gt lt sample in k space gt lt lattice gt 100 lt lattice gt GJS
180. mensions These are specified with the lt moment_group gt element We do not have one of these in our worked example but as an example we might need the total number of atoms lt moment_group gt lt moments gt number lt moments gt lt integrate_dimension gt yes lt integrate_dimension gt lt CDATA number phixphi gt lt moment_group gt The child elements here name the moments to be defined and then describe which if any of the transverse dimensions are to be integrated The CDATA block defines the moments them selves The number moment in this example will be integrated over the only transverse dimension lt moment_group gt and lt functions gt elements may be used any number of times in any order The propagation of the integration equations themselves is performed where the lt vectors gt tag is placed Something else to keep in mind is that when a smooth function is required it is often tempting to use transcendental functions such as sqrt exp sin cos log However these functions are computationally expensive particularly on processors that were not de signed with such functions in mind If they are only used once to pre calculate a constant vector as is done here then fine but if you include them in your main equations then ex pect a big reduction in performance If a smooth function that depends on the propagation dimension is required then unless it is essential to the model it is usua
181. ment with which one can set the output precision to either single or double precision The syntax for this is as follows lt output format binary precision double single gt SL mozo nds rase gt lt output gt where double single means the options are either double or single with double being the default option Notice that the format attribute is also set to binary this is to emphasise that it is pointless specifying the precision without the format since the precision attribute is meaningless for ascii output Using this option and rerunning the atomlaser simulation we find that the file size of atomlasermg0 dat is 168 kB and atomlaser xsil is 4 kB which overall is 21 of the size of the original ascii output 2 6 Initialisation of field vectors from file In Chapter 1 Section we initialised the field inside the lt vector gt element by using C C code It is also possible to already have these vectors calculated and stored in a file which xmds can then load and use to initialise the field This feature can be useful if the calculation of the vectors is particularly difficult and you don t wish for xmds to have to calculate them or you may have already generated the data from another program and so going through the hassle of getting xmds to recalculate the data is a waste of time Anyway it can be handy to do on some occasions and so xmds provides a means for you to do this via 2 6 INITIALISATION OF FIELD V
182. n 1 above provided that you also meet all of these conditions a You must cause the modified files to carry prominent notices stating that you changed the files and the date of any change b You must cause any work that you distribute or publish that in whole or in part contains or is derived from the Program or any part thereof to be licensed as a whole at no charge to all third parties under the terms of this License c If the modified program normally reads commands interactively when run you must cause it when started running for such interactive use in the most ordinary way to print or display an announcement including an appropriate copyright notice and a notice that there is no warranty or else saying that you provide a warranty and that users may redistribute the program under these conditions and telling the user how to view a copy 233 of this License Exception if the Program itself is interactive but does not normally print such an announcement your work based on the Program is not required to print an announcement These requirements apply to the modified work as a whole If identifiable sections of that work are not derived from the Program and can be reasonably considered independent and separate works in themselves then this License and its terms do not apply to those sections when you distribute them as separate works But when you distribute the same sections as part of a whole which is a work based on the Pr
183. n be confident that discretisation error is not a significant problem in the simulation output Once you are sure that your simulation is behaving nicely you can then turn off error checking by setting error check to no thereby speeding up the simulation This is especially important for those whose simulations are going to take a long time to run So test your simulation to make sure that the error is low for a short simulation run and then for the main run turn error checking off 2 2 MPI automatic parallelisation of simulations One of the most powerful features of xmds is its ability to automatically parallelise sim ulations We go through an extended discussion of this with specific focus on stochastic simulations in Chapter 4 but it is of worth mention here as well Not only can xmds paral lelise stochastic simulations by running each stochastic path on a separate computer but it is also able to parallelise the computation of deterministic problems as well It does this with the help of a package known as MPI which stands for the Message Passing Interface and is a means to organise the communication and computation associated with parallel simulations To parallelise your simulation all you have to do is add one line of code Honestly We re not joking To turn this feature on in your code you need to add the line use mpi yes use mpi and that s it What xmds will do when running your simulation is split up the computati
184. nctions available is the same as that for C programming as well as some complex functions defined in the xmdscomplex h header file in the source directory The ones most likely to be of use are summarised in Table 8 4 And remember avoid using transcendental functions in the main integration equations if speed is important 134 FUNCTIONALITY Example for lt prop_dim gt z lt prop_dim gt lt noises gt 2 lt noises gt lt dimensions gt x y lt dimensions gt Variable Available in Z lt integrate gt lt k_operators gt lt filter gt lt sampling gt zor kz lt post_propagation gt dz lt integrate gt lt k_operators gt x or kx lt vector gt lt integrate gt lt k_operators gt yorky lt sampling gt dx or dkx and lt post_propagation gt where dimension dy or dky hasn t been collapsed ni n2 lt vector gt integrate filter Table 8 5 Automatically declared variables It is also possible to create piece wise functions by using if statements which must be done in the usual C syntax Here are two alternate syntaxes that will turn on functions dampingi and damping2 after t 1 0 and t 2 0 respectively double damping if t gt 1 0 dampingi 1 0 ellise dampingi 0 const double damping2 t gt 2 07 1 0 xmds automatically declares a number of variables depending on the element for the equations to reference These are summarised in Table If a
185. nd several hours of computer time for the simulation to come back to you and say I m done and you haven t got any results Fortunately xmds doesn t let you write a simulation without actually specifying any output Therefore we need to tell xmds what output we want from the simulation and to do this we use the lt output gt element The lt output gt element is just a container for the other tags that specify what is to be output The two tags that are contained within the lt output gt element are lt filename gt and lt group gt 14 STARTING FROM SCRATCH The lt filename gt tag fairly obviously specifies the filename of the output data file This tag is optional and defaults to the simulation name i e the value of lt name gt directly within the lt simulation gt tag with the string xsil appended For example if we didn t specify a filename for the simulation we ve created here then the filename xmds would use would be basicSim xsil The output data file is in the XSIL format 2 which is a handy interchange format also using XML The lt group gt tag contains a description and to a degree the definition of the moments of the output data which can be such things as the power density of the field s or just means and standard errors of the field s More than one output group can be specified but at least one must be given Post processing may be performed before an output group is written to file allow
186. ne of code Part II Numerical Modelling Theory TT Introduction Unfortunately appropriate mathematical models for the majority of natural phenomenon generate equations that are non integrable This means that the solutions to the problem cannot be written as exact analytic expressions and those who would seek solutions must appeal to the modern computer to calculate them by a numerical procedure The capability of modern computers continues to grow at a phenomenal rate and hence so also does the field of computational physics In fact it has been described by some as a third way of doing physics in addition to the more standard theoretical and experimental approaches Compu tational physics represents an unusual marriage of these two approaches which traditionally have always been rivals in that the governing equations must be based on believable theory yet each simulation is very much like an experiment to determine what actually happens However due perhaps to its rapid growth computational physics is in need of improvement in at least one particular area In a recent review titled Microscopic simulations in physics 4 Ceperley writes Sadly the lore of experimental and theoretical physics has not yet fully penetrated into computational physics Before the field can advance certain standards which are common place in other technical areas need to be adopted so that people and codes can work together Here Ceperley is add
187. nitial const double zo 20 0 initial alee lt globals gt lt field gt lt name gt main lt name gt lt samples gt 1 lt samples gt lt vector gt lt name gt main lt name gt lt type gt double lt type gt lt components gt x y z lt components gt lt CDATA x ES xo yo by Alejandro L Adapted from the example Garcia conditions conditions conditions 16 STARTING FROM SCRATCH gt lt vector gt lt field gt lt sequence gt lt integrate gt lt algorithm gt RK4EX lt algorithm gt lt interval gt 10 lt interval gt lt lattice gt 10000 lt lattice gt lt samples gt 200 lt samples gt lt CDATA dx_dt sigma y x Gl ux c y URS dz_dt x y b z 11 lt integrate gt lt sequence gt lt output gt lt sampling gt moments x0ut yOut z0ut lt moments gt lt CDATA xDut x yOut y zu 225 Ie lt sampling gt lt output gt lt simulation gt Here is a link to the finished gzipped script file on the xmds web site http www xmds org 1 1 7 Making the simulation and getting results Now that the simulation script is ready it is just a matter of getting xmds to generate the C source code and compiling that with your system s C compiler This is a very simple process and in the vast majority of cases all one has to do is enter the following at the command prompt xmds lorenz xmds A file called l
188. ns gt lt integrate gt lt sequence gt lt simulation gt 11 22 2 6 min time step optional lt min_time_step gt double lt min_time_step gt Contains double Subelement of Path to tag lt simulation gt lt sequence gt lt integrate gt lt min_time_step gt Description The minimum time step an integration algorithm should try before halting that pass of the integration This option applies to the adaptive algorithms ARK45 and ARK89 only The default value is 1e 13 that is 1071 A value of O disables the check on the time step which will slightly speed up integrations that do not require the check Example lt simulation gt lt sequence gt lt integrate gt lt min_time_step gt 1e 10 lt min_time_step gt 202 LANGUAGE REFERENCE lt integrate gt lt sequence gt lt simulation gt 11 22 2 7 smallmemory optional lt smallmemory gt bool lt smallmemory gt Contains boolean Subelement of Path to tag lt simulation gt lt sequence gt lt integrate gt lt smallmemory gt Description Defaults to no In this case when using the ARK45IP algorithm xmds calculates and stores six copies of the derivative operators per step as each of them can be reused once For problems with memory limitations it might be better to assign yes to this element then each array of derivatives needs to be calculated twice per step but significantly less memory should be used E
189. ns gt tags since they are all required just for the transverse dimensions We want to sample the first point so we set lt samples gt to 1 Our variable z is complex so we leave the first two tags of the lt vector gt block as is and set the lt components gt tag to z which we don t define in Fourier space and we just set to the initial value zo This gives a lt field gt block of GIE Field to be integrated over gt lt field gt lt name gt main lt name gt lt samples gt 1 lt samples gt ASS pS TP OCIO TN gt lt vector gt lt name gt main lt name gt INotice that we don t actually use the word gamma here This is because gamma is often used in maths libraries for such things as the gamma function etc so we try and avoid name conflicts as much as possible and use the word gam 2we could have just set z to some number later in the setup of the field however this way allows us to change the initial value of z later without having to dig through too much code and we can easily see how to replace the lt globals gt block with an lt argv gt block and so have more dynamic control over the variables being put into the simulation 4 1 WITHOUT MPI 67 lt type gt complex lt type gt si cama uy 0 v e c toi gt lt components gt z lt components gt lt names of components gt lt fourier_space gt no lt fourier_space gt lt defined in k space gt lt CDATA ZO 1 gt lt vector g
190. ntains a tabulated summary of the syntax 9 1 The simulation element lt simulation gt lt name gt nlse lt name gt 138 WORKED EXAMPLE NLSE XMDS lol N CO OC CUA IN 09 NOS N ju 9 WAWA SE FIGURE 9 1 Results for nlse xmds lt prop_dim gt z lt prop_dim gt lt error_check gt yes lt error_check gt lt stochastic gt no lt stochastic gt SI MO SS gt lt simulation gt The root XML element of the input file must be a lt simulation gt element which contains the high level description of the problem Within the lt simulation gt element xmds will first look for a simulation lt name gt This information is optional the default being the file name minus its last extension It is used as the base for the output cc file and is the name of the executable that xmds will make The lt prop_dim gt assignment defines the name of the main propagation dimension It is compulsory as it is expected that whichever symbol is used here will also appear in the equations the user supplies The lt error_check gt assignment is optional the default being yes This option tells xmds whether to write the code so as to repeat the simulation with half the step sizes in the integration segments and subsequently compare the two sets of results It is a good idea to leave this option on until confident that the errors are acceptable for the simulation in question
191. ntegrate element lt integrate gt lt algorithm gt RK4IP lt algorithm gt lt interval gt 20 lt interval gt lt lattice gt 1000 lt lattice gt lt samples gt 50 lt samples gt lt k_operators gt lt constant gt yes lt constant gt lt operator_names gt L lt operator_names gt lt CDATA L rcomplex 0 ktx kt 2 1 1 gt lt k_operators gt lt vectors gt main vcl lt vectors gt lt CDATA dphi dz L phi i phi phi phi phi damping In lt integrate gt Here lies the heart of what xmds is all about and not surprisingly it forms the most com plex element The lt algorithm gt assignment is optional and will default to SIEX for stochas tic simulations and RK4EX for non stochastic simulations The assignments lt interval gt lt lattice gt and lt samples gt are all compulsory and respectively represent the integration range total number of steps and number of samples for each output moment group to take within these steps Thus each integer in the lt samples gt assignment must be either zero or else a factor of lt lattice gt and if not xmds will exit with an error message to that effect The rest of the lt integrate gt element consists of writing down the form of Equa tion in a high level format Here the algorithm is RKAIP which as explained earlier is an interaction picture algorithm The dispersion term in Equation is represented here by the linear operator L acting on field c
192. ntegration Now suppose the problem uses N noises and has a transverse lattice of M points Within the SIEX RK4IP and the RK4EX integration algorithms the main field vector is swept through severals times in the course of each time step thus a M x N vector of noise must be calculated at the beginning of the time step and referenced during the calculation of the main field vector s derivatives In the full step case two such vectors are calculated and 10 3 FIBRE XMDS 155 then averaged to provide the equivalent full step noise However the SIIP algorithm only sweeps though the main field vector once for each time step and so only N noises need be calculated at a time thus saving on memory and RAM access provided two independent random number generators are used the first generator is used for the first half step the second for the second half step and the average of both is used for the full step In C this is done with the erand48 n function which uses the 48 bit integer n to generate the next random number advancing n to the next in sequence in the process The states of the two independent generators are simply two independent integers n and nz The user supplies the initial values for these integers in the lt seed gt assignment Stochastic problems are very well suited to parallel computer architectures as different paths can run on different processors and do not have to transfer information until the integration is complete
193. nto account contributions from higher order derivatives In practice these methods are used rarely due to the availability of a method based on the Fourier Transform The partial derivate with respect to dimension 2 is calculated quite simply using Equation 6 32 O o 1 ik x aa Sai oat ee 1 dke ik k a J kiya k 6 32 This method can equally be applied to discrete data Firstly a discrete Fourier Transform is applied to the field data az then the data is multiplied point by point by the corresponding 6 4 NUMERICAL METHODS FOR PARTIAL DIFFERENTIAL EQUATIONS 89 ik values and then the discrete inverse Fourier Transform is applied to obtain the a data Extending this method for higher order derivatives is trivial Equation 6 33 illustrates some example mappings J k ku JE f k k 2 ox 2 OxOy ika k gt ike iky k k 6 33 This method becomes computationally more efficient than the matrix method when a number of higher order derivatives are desired for example when a Laplacian operator is being applied in multiple transverse dimensions The main limitations of this method are that it implies periodic boundary conditions and that it cannot be implemented on a non regular grid 6 4 2 Explicit Picture Methods In the explicit picture all linear and nonlinear operators are explicitly calculated and summed to generate the total derivative for each field
194. ogram the distribution of the whole must be on the terms of this License whose permissions for other licensees extend to the entire whole and thus to each and every part regardless of who wrote it Thus it is not the intent of this section to claim rights or contest your rights to work written entirely by you rather the intent is to exercise the right to control the distribution of derivative or collective works based on the Program In addition mere aggregation of another work not based on the Program with the Pro gram or with a work based on the Program on a volume of a storage or distribution medium does not bring the other work under the scope of this License 3 You may copy and distribute the Program or a work based on it under Section 2 in object code or executable form under the terms of Sections 1 and 2 above provided that you also do one of the following a Accompany it with the complete corresponding machine readable source code which must be distributed under the terms of Sections 1 and 2 above on a medium customarily used for software interchange or b Accompany it with a written offer valid for at least three years to give any third party for a charge no more than your cost of physically performing source distribution a complete machine readable copy of the corresponding source code to be distributed under the terms of Sections 1 and 2 above on a medium customarily used for software interchange or c Accompany it
195. oment group filter 240 004 medo nomo ae de 210 11 22 3 2 functions filter Leu ack ee des eae oe que eh cee ud 210 11 22 3 3 vectors filter ue ded eee a 211 11 22 3 4 fourier space filter ocu wa we oe we R x TR Ree 211 uns sss s eae need 212 11 22 3 5 1 prop_dim cross propagation 4 44 4 44 40 212 11 22 3 5 2 vectors cross propagation 213 DL223 4 OS a Xo wed 4e uod ox a SNS EEE BS OR 4 213 CONTENTS 11 22 4 1 filename breakpoint 11 22 4 2 fourier space breakpoint 11 22 4 3 vectors breakpoint 11 23 2 1 1 type sampling 11 23 2 1 2 fourier_space sampling 11 23 2 1 3 vectors sampling 11 23 2 1 4 lattice sampling 11 23 2 1 5 moments sampling 11 23 2 2 1 fourier_space post_propagation 11 23 2 2 2 moments post_propagation A The XSIL input output syntax B The xsil2graphics utility program lt gt y 5 O 2 x C loadxsil m utility script RS bate A CR nd pit Rea i occ AA d sub RP med AS rrr MARE uoo eben koe Od c Qux need DP oh ne de Dp up e RE die rios eie Seius de Sah Bose bue a RW ERROR do zoe ca os a dex ens DU Da we do as LES Copyright seora box dex vm bom deo hok drew E Ra dk d Read D Gnu General Public License VI Bibliography and Index Bibliography LIQUID 2 s s De v wow a e oe a Boke et ROOD ROR e W AS 11 23 1 filename output J o E ee EE Reds bie Ge a edo ee de Pn de do D i 11 23 2 1 sampling 4 4
196. omponent phi The presence of the term L phi in the equations causes the field to be evolved in Fourier space in accordance with the algorithm listed in Section The W operators are everything left over once terms such as this are removed 9 6 THE INTEGRATE ELEMENT 143 The lt k_operators gt element contains just what its name implies k space operators Within this element a lt constant gt declaration is used to indicate whether or not these op erators depend on the propagation dimension This declaration is optional with the default being no but this will produce slower code If the operators do not depend on the propa gation dimension then you will get faster code is you enter lt constant gt yes lt constant gt An lt operator_names gt assignment is compulsory One or more operator names may be listed here and they all must be explicitly defined in the code Refer to Equation for how to write the LY operators as functions of the Fourier space dimensions The operator names must be valid names as far as the C language is concerned sce Section 8 3 If there are no non zero linear operators then the entire lt k_operators gt element may be omitted Because of the way xmds works internally it is possible when using an interaction picture routine to write the main equations in such as way as that everything works fine but the results will be incorrect The problem is that xmds does not closely inspect the main
197. on 1 0 gt lt simulation gt lt name gt kubo_tutorial lt name gt lt the name of the sim gt lt author gt Paul Cochrane lt author gt lt the author of the sim gt lt description gt Kubo oscillator example simulation This formally represents a simple model of an oscillator with a mean frequency omega perturbed by a noise term xi t Adapted from the example given in Handbook of Stochastic Methods C W Gardiner 1997 lt description gt lt Global system parameters and functionality gt lt prop_dim gt t lt prop_dim gt S name Gr manim paso papi won im gt lt stochastic gt yes lt stochastic gt lt defaults to no gt AE these three tags only necessary when stochastic is yes gt lt paths gt 1024 lt paths gt Gl mo of parns gt lt seed gt 1 2 lt seed gt ASES E SS O E nc we gen lt noises gt 1 lt noises gt Wc no of noises gt lt use_mpi gt no lt use_mpi gt lt de F am S vo mo gt 4 1 WITHOUT MPI 69 lt error_check gt yes lt error_check gt amp l dqeta lts to yes gt lt use_wisdom gt yes lt use_wisdom gt lt i defaults to ug gt lt benchmark gt yes lt benchmark gt SI defaults to mo gt lt use_prefs gt yes lt use_prefs gt lt I gt detavdatis LO yes lt l Global variables fon the simulation gt gt lt globals gt lt CDATA double omega 0 mean frequency the word gamma is already used in
198. on as in Chapter 1 Section 1 1 3 Although xmds will continue initialisation even if it cannot find all the variables in a vector in the XSIL file it will not continue if it cannot find any variables it will print a warning about any variables that it cannot find in the XSIL file Note that the sequence of initialisation steps for each element in a vector is to first initialise the element to zero then to use any code in the CDATA section if present and finally to initialise from the XSIL file if the variable is present in the moment group Hence initialisation from the XSIL file will override any initialisation in the CDATA section The dimensions of a vector can be initialised in any combination of z space and k space by using the lt fourier_space gt tag in the same way as it is used for initialisation from C C code however the default is that each dimension is initialised in x space There are some restrictions on the geometry of the moment group in the XSIL file however these conditions depend on whether the geometry matching mode specified by the attribute geometry matching mode of the lt filename gt tag is set to strict mode or loose mode In strict mode the following conditions apply 1 The moment group must have the same number of dimensions as the field In other words the moment group can only have been sampled once as sampling a moment group a number of times introduces an extra dimension the propagation dimension
199. on of the field and pass these parts of the overall computation to different computers where it is solved faster than possible by doing so on a single machine One good reason for splitting a simulation like this up and processing each part of the field on different processors is because for some very large simulations the memory requirements are too large and therefore won t fit on one computer using MPI is the only way to solve the problem There is one major caveat here however ONLY use MPI for deterministic problems on a supercomputer or a cluster setup where there is very small network latency it doesn t take long for computers to talk to one another This is very important because the Fourier transforms require a lot of communication and if the network between nodes of the cluster is slow then this will reduce the speed of the computation significantly probably making it faster to run on a single cpu Nevertheless if you do have access to powerful computing facilities then by all means use this feature For stochastic problems there is a second option you may wish to use which changes the way the different paths are allocated between processors This can be altered by using the MPI Method tag which can take the values Scheduling or Uniform 2 3 BENCHMARKING 35 2 3 Benchmarking To get an idea as to how long your code is going to run when you scale the various simulation parameters up for a long simulation or ju
200. one vector must be called main Example 194 LANGUAGE REFERENCE lt simulation gt lt field gt lt vector gt lt name gt main lt name gt lt vector gt lt field gt lt simulation gt 11 21 6 2 filename vector optional lt filename gt string lt filename gt Contains string Attributes optional format ascii binary xsil Attributes for the XSIL format optional moment_group N optional geometry matching mode strict loose Subelement of Path to tag simulation field lt vector gt filename Description Tells xmds the file from which to load the initial field The lt filename gt tag accepts one optional attribute format This attribute can take one of three options ascii binary or xsil where ascii is the default option However xsil is the most robust of these formats as it can load any binary XSIL file produced by XMDS irrespective of what architecture the file was produced on See Chapter 2 Section 2 6 1 for more information on loading XSIL files Caution should be taken here when using the binary format since loading a binary file is more difficult to get correct than loading an ascii file However the added difficulty is offset by being able to have smaller and more complex input data A first thing to note is that the byte ordering of the system you are going to be loading the file into has to be the same as the file itself For instance creating
201. onvergence of a numerical method will be less than or equal to half the order of the deterministic convergence 6 3 2 The Improved Euler Method Here two evaluations of the derivative f t a are made The first one is based on the current value a and is used to estimate the next value of a as in the explicit method This estimated next value is then used to calculate the next value of the derivative which in turn is used to make a second estimate of the next value of a as in the implicit method These two estimates are then averaged a a h fl an 6 15 a a th flint a 6 16 1 Anti 5 an aia 6 17 6 3 NUMERICAL METHODS FOR DIFFERENTIAL EQUATIONS 85 The two estimates of f t a account for the second order derivative of a so the error per step is reduced to order h and the integral is correct to order h This method is also known as the improved Euler Cauchy method or Huen s method It is also known as the second order Runge Kutta method Unfortunately this method is not guaranteed to converge for stochastic integration when the noise terms are naively added to the function on the RHS 6 3 3 The Semi Implicit Method This method is very similar to the improved Euler method Essentially the method lies half way between the explicit and implicit methods in that the new value for a is calculated based on the derivative at the mid point of the step Similarly to the implicit method the current value of a
202. ore complex than is yet possible with xmds The first thing xmds does is run the input file through its own XML parser If the input file contains bad XML syntax then these will be the first errors to be picked upon xmds then processes the input file and writes the code dedicated to solving the particular simulation and compiles it xmds nlse xmds compiling gcc pthread 03 ffast math funroll all loops fomit frame pointer o nlse nlse cc I home cochrane bin lstdc 1m lxmds L home cochrane bin 1fftw_threads 1fftw nlse ready to execute All that remains is to execute the compiled program nlse Refer to the kubo xmds example in Section for execution of a parallel problem Or see Chapter 8 1 3 Preferences As of xmds 1 3 1 it is possible for people to specify the simulation build options in a preferences file Previously if one wanted to change how the simulations were built one had to recompile xmds with the relevant compile flags etc Now all one has to do is specify the compilation flags in a preferences file xmds looks for a file called xmds prefs in HOME xmds and if not found there in the directory local to the xmds script you are trying to compile The format of the preferences file is lt compile flag name gt lt compile flag value gt For example XMDS_CC gcc Comments can be added by using a hash character Any thing after and includ ing the hash are ignored As another example of a prefe
203. orenz should now appear in the same directory as your simulation script this file is the simulation binary executable file To run the simulation merely execute the binary file by entering its name at the command line You should see something like this lorenz Beginning full step integration 1 1 A BASIC SIMULATION 17 Sampled field for moment group 1 at t 0 000000e 00 Sampled field for moment group 1 at t 5 000000e 02 lt snip gt Sampled field for moment group 1 at t 9 950000e 00 Sampled field for moment group 1 at t FOO 000 Oe On maximum step error in moment group 1 was 1 248578e 05 Once the program has finished running you should find in the same directory as the binary executable a file called lorenz xsil This is the file containing your output data in a handy XML based format that can be used to interchange data between various other formats We ll look at the output here in two programs namely Matlab or Octave and Scilab Matlab is a commercial numerical programming language and environment which has very powerful graphics capabilities and is used in the scientific community extensively Octave is a free program which is highly compatible with Matlab Scilab is very similar to Matlab however it is free to download and install but doesn t have quite the same quality as that made by Matlab Nevertheless Scilab is free and is a handy alternative if your budget can t stretch to Matlab There are subtle differ
204. ors initialised from file are assumed to be entirely in normal space 9 5 The sequence element lt sequence gt lt integrate gt lt More xmds code gt lt integrate gt lt More xmds code gt lt sequence gt Next comes the lt sequence gt element which can be used in two different contexts The first context as used here is the top level lt sequence gt element and is simply for defining the sequence of segments that happen to the field after initialisation In a stochastic simula tion this sequence is repeated for each integration path but in this non stochastic example the sequence is only stepped through once In this context xmds expects to find any number 142 WORKED EXAMPLE NLSE XMDS of integrate lt filter gt or sequence elements which leads us to the second context In the second context lt sequence gt elements may be used as child elements within parent lt sequence gt elements in which case they represent a sub loop to be repeated one or more times The number of repetitions is defined by a simple lt cycles gt assignment within which is expected to be the first element In the absence of a lt cycles gt element the number of cycles defaults to one but if it is present it must come first Any number of lt integrate gt lt filter gt or lt sequence gt elements may then follow Obviously the ordering of segments within any sequence element is important 9 6 The i
205. ound at Version 2 1 x and optionally 3 x of the fftw library is required for proper execution of xmds xmds expects that the fftw h header file is in one of the usr include or the usr local include directories with the libraries in a sibling usr lib or usr local lib directory These are the standard locations If the fftw installation is non standard then you will need to use the with fftw path option for the configure program pointing to the parent directory containing the fftw include and lib sub directories If parallel processing is desired FFTW must be built with MPI or OpenMP options enabled The main method of installation is from a tarball Download xmds 1 6 5 tar gz and 8 1 INSTALLING AND RUNNING xmds 127 then as root tar xzvf xmds 1 6 5 tar gz cd xmds 1 6 5 configure make make install or as a user to be installed in the bin directory of your home directory tar xvzf xmds 1 6 5 tar gz cd xmds 1 6 5 configure with user make make install For help on the various configuration options one can run the command configure help which will show a very long list of options which can be changed to customise the way xmds is installed and the options that it can use to build simulations The different variables that one can change are e XMDS_CC the C compiler used by xmds to compile simulations e XMDS_CFLAGS the C compiler flags used to build simulations e XMDS_LIBS the libraries
206. our NIS z t zoe 2 10 10 Since the variance of the noise terms scale with the inverse of the integration step size the integration method suffers a loss of order with regard to error vs step size While this is normally a second order method for non stochastic problems it becomes a first order method for problems with noise as was explained in Section 6 3 1 10 3 fibre xmds lt xml version 1 0 gt Si emae aloe OLSE See iom lt simulation gt real z 0 5 0 5 156 MORE EXAMPLES Eg lt real z gt 2 4 6 8 10 t a Single path 0 8 0 2 0 b 1024 path mean FIGURE 10 2 Results for kubo xmds lt name gt fibre lt name gt lt prop_dim gt t lt prop_dim gt lt error_check gt yes lt error_check gt lt stochastic gt yes lt stochastic gt lt use_mpi gt no lt use_mpi gt lt paths gt 1 lt paths gt lt seed gt 1 2 lt seed gt lt noises gt 2 lt noises gt lt globals gt lt CDATA const double ggamma 1 const double beta sqrt 2 2 M_PI ggamma 10 1 gt lt globals gt lt field gt lt name gt main lt name gt lt dimensions gt ps lt dimensions gt lt lattice gt 50 lt lattice gt domains 5 5 lt domains gt lt samples gt 1 lt samples gt lt vector gt lt name gt main lt name gt lt type gt complex lt type gt lt components gt phi lt components gt 10 10 3 FIBRE XMDS 157 lt fourier
207. output gt lt noises gt lt paths gt lt benchmark gt obsolete use visdom gt obsolete kuse prefs gt Subelement of None Path to tag lt simulation gt Description Container tag for the xmds simulation code Example lt simulation gt lt U senes CERS gt lt simulation gt 11 2 name simulation optional lt name gt string lt name gt Contains string 178 LANGUAGE REFERENCE Subelement of Path to tag simulation lt name gt Description The name of the xmds simulation Defines the name of the output C code file and subsequently the name of the simulation binary executable This tag is optional and if not specified then the output filenames will be derived from the xmds script filename minus its last extension For instance for the atomlaser xmds script the base for the output cc file will be atomlaser This isn t really a very instructive example Example lt simulation gt lt name gt atomlaser lt name gt lt simulation gt 11 3 prop dim simulation required lt prop_dim gt string variableName lt prop_dim gt Contains string Subelement of Path to tag lt simulation gt lt prop_dim gt Description The name of the main propagation direction This name must appear in the equations supplied This condition however is not checked and therefore a possible source of error for the user Example lt simulation gt lt prop_dim gt z lt prop
208. ow we come to the sequence section of the code There s only one such sequence in this simulation so we really only need to worry about the lt integrate gt block The algorithm we ll use here is the fourth order Runge Kutta method in the explicit picture We want the wave to circle the system the once so we ll choose the lt interval gt and lt lattice gt so that it does so while having a step size of about 0 002 We therefore set lt interval gt to 1 and lt lattice gt to 500 Only 50 points are really necessary for the output of this hence lt samples gt is set to 50 The next section to define is the lt k_operators gt element The k space operator is con stant in time and we ll call it L which is sort of a convention in the xmds community Since the operator acts in Fourier space the spatial derivative merely becomes a multiplication and so we can set the CDATA block within the lt k_operators gt element to lt CDATA L rcomplex 0 0 v kx gt Note that we have used the rcomplex function and that the operator is complex This is because of the mapping of derivatives in z space to k vectors in k space remember that 2 ee 3 5 There are no lt vectors gt to define so we can just delete this line however we must write down the differential equation we re trying to solve So remembering Equation 3 1 and the fact that xmds uses a special syntax for writing the differential equation in the lt integr
209. pend limf to XMDS_LIBS this adds the icc native support for its maths libraries Debugging By default xmds is configured to use quite aggressive optimisations when compiling simulations If however you suspect something is going wrong and you wish to debug the simulation binary directly using gdb dbx or another symbolic debugger then you will need to put symbolic debugging information into the binary executable To do this replace the default options by setting the XMDS_CFLAGS variable to g Profiling There may be instances when one wants to find out what part of the code is taking up the most time when running a simulation This is generally speaking a part of debugging and testing a simulation and not normally part of using xmds However if you re interested in seeing what lines of code are using the most time you ll want to add profiling information to the code To do this you will need to add either the p or gp option to the XMDS_CFLAGS variable The p option generates extra code for profiling with the prof utility and the gp option generates extra code for profiling with the gprof utility 50 EXTRA AND ADVANCED FEATURES 2 9 Breakpoints Breakpoint elements are parts of a simulation similar to an lt integrate gt or a lt filter gt element that cause the state of some vectors of the simulation to be saved to an XSIL file when the breakpoint element is hit This can be used early in a simulation to enable you to c
210. phy and Index 237 Bibliography A L Garcia Numerical methods for physics Prentice Hall 1994 XSIL documentation URL http www cacr caltech edu SDA xsil C W Gardiner Handbook of Stochastic Methods Springer Verlag Berlin 1997 2nd ed D M Ceperley Rev Mod Phys 71 438 1999 E Kreyszig Advanced Engineering Mathematics John Wiley amp Sons 1999 P D Drummond Comp Phys Comm 29 211 1983 B M Caradoc Davies URL http www physics otago ac nz research uca resources bmcdthesis html M J Davis URL XML syntax specification URL http www w3 org TR 2000 REC xm1 20001006 E R Harold and W S Means XML in a Nutshell O Reilly amp Associates Inc Se bastopol USA 2002 2nd ed Document Object Model specification URL http www w3 org DOM 239
211. pled in normal space xmds uses a coarser lattice that still spans the entire domain whereas if it is in Fourier space it samples a narrower window of the k space domain centered on the zero momentum point k 0 You may of course specify the same number of lattice points here as there are in this dimension in which case the output moments are simply mapped point for point If there are no transverse dimensions xmds will not look for either of the above two assignments Next comes the lt moments gt assignment which lists the names of the moments Any number of moments are possible provided there is at least one Finally the code used for calculating these moments is inserted as the content of the lt sampling gt element The moments are complex variables and so the code may be written accordingly Next a lt post_propagation gt element may optionally follow the lt sampling gt element In this element there must be a lt fourier_space gt element with a yes or no entry for the propagation dimension and as many more as there are remaining transverse dimensions in the sampled output In the case of partially collapsed output the correspondence is sequential For example if the field had three transverse dimensions x y z and the y dimension was collapsed by using a lt lattice gt assignment within the lt sampling gt element like 50 1 50 then a lt fourier_space gt assignment of yes no yes would mean transfor
212. plying C code Here we are doing the latter In the absence of lt filename gt assignment xmds looks for a lt fourier_space gt assignment followed by the initialisation code The lt fourier_space gt assignment must contain as many yes or no entries as there are transverse dimensions in the field which in this case is only one A yes entry means the corresponding dimension is in Fourier space when initialised and vice versa When initialising from code the code is included as the content of the lt vector gt element There must be an equation for each component and these may include functions of lt global gt variables and of the transverse dimensions Again see Section for more information Two points to remember firstly if the field components are complex variables they may be initialised using the complex functions in Table and secondly every field component should be explicitly initialised even if it is only being initialised to zero Finally when initialising from file xmds expects that the file is simply a sequence of numbers in ASCII format with as many numbers as the product of lattice points and vector components In usual C convention the component index varies most rapidly then the right most transverse dimension lattice and so on The left most transverse dimension lattice dimension varies most slowly If the vector is complex then each data point is read in as a sequential real and imaginary pair All vect
213. ptional lt noises gt int lt noises gt required if lt stochastic gt is yes Contains integer Attributes optional kind gaussian gaussFast poissonian uniform mean int re quired for poissonian Subelement of Path to tag lt simulation gt lt noises gt Description The number of noise terms in the simulation xmds automatically declares the variables n_1 n_m where m is the number of noises specified As of xmds 1 3 3 the lt noises gt tag now accepts attributes Currently these are gaussian gaussFast poissonian and uniform Each refers to a different noise source gaussian refers to the original noise routine that xmds uses It produces a Gaussian distributed variable with mean of zero and variance of 1 product of the propa gation and transverse step sizes gaussFast is a slightly faster implementation of the same routine and will eventually replace the old gaussian routine 182 LANGUAGE REFERENCE poissonian is a Poissonian noise source and requires the mean_rate attribute to be set to the desired mean rate of the Poissonian noise distribution Simulations where the mean rate is a function of the propagation dimension or other variables are not supported via this method but judicious use of the functions element can allow this rate to be adjusted manually uniform gives uniformly distributed noise on the interval 0 1 Example lt simulation gt lt noises kind gaussian gt 1 lt noise
214. r in this case y then the first element of the third vector z here and then the second element of the first vector and so on One way of describing this is in terms of C C code The data is expected in this format x 0 y LO z 0 x 1 y 1 z 1 and so on until the end of the data Another way of describing this is in terms of the actual data and so here is how the file input dat will look 2 0 5e 2 10 1 0 1e 3 20 0 0 1e 5 30 1 0 2e 4 50 2 0 Te 4 1e3 If this seems unnecessarily complicated it is However this is the way the data is expected and so we have to behave the way xmds expects otherwise our simulation will not work properly As it turns out storing the data within memory in this fashion means that calcu lations are performed on contiguous blocks of memory and therefore are a lot faster than if entire vectors were stored with their elements next to one another This is a significant point for the memory utilisation internal to the simulation and for maintaining the speed of xmds 2 6 INITIALISATION OF FIELD VECTORS FROM FILE 43 simulations However it may be possible in future versions of xmds for the input data to be specified more logically i e have x defined first then y etc and then for the simulation to reorganise the data internally so that calculations are performed efficiently and quickly If your input data is binary instead of ascii as we have above then you would use the
215. r nonlinear partial derivatives that may 123 124 FUNCTIONALITY have spatial dependence for example L 23 Equation is Schr dinger like in that it is first order in the propagation dimen sion As an example if the field was a vector field representing three dimensional fluid flow then there would be three components to the field which itself would exist in four dimen sional space with perhaps the propagation dimension being time and the three transverse dimensions being space xmds can only handle problems whose only non local terms are integrals over some or all of the transverse dimensions For local equations the evolution of the field at any point in x is only ever a function including derivatives of the field components at that exact point in space not anywhere else This is clearly restrictive of the range of physical problems that may be modelled However many non local problems are only so because of variables that propagate in a dimension which is transverse to the main propagation dimension and thus are easily modelled by including a secondary field which is propagated orthogonally along this transverse dimension This is the purpose of the b field which is the secondary or cross field propagating in the transverse dimension z where c Z 0 In order to maximise the efficiency of the generated code while being able to handle a broad range of problems xmds utilises a range of algorithms These were covered in Sec
216. r operators and the b field is evolved over the time step using the nonlinear operators and the semi implicit algorithm Thus the algorithm becomes quite simply Sh z9 ki l a e a 2 Calculate if required 3 39 x h 4 For each point in space do a arsa b For N 1 iterations do i a ar iA N x a c a 2 ar h N x a 5 39 x h IhC a9 k 6 a e2 a The first and the last steps are easily performed with the field in Fourier space since the linear operators reduce to k multipliers Note that the dependence of the linear operators L x k i on the propagation dimension ought to be weak in comparison to the size of the time step If they do vary with time then the exponentials will have to be calculated at every time step which can be computationally expensive with some CPU architectures Otherwise the exponentials may be pre calculated and tabulated for reference Further advantages of this method are that the forward and backward Fourier transforms may be performed in place i e within the current memory block as can the point by point multiplication of the field by the exponentials Thus no extra copies of the field are required 6 4 NUMERICAL METHODS FOR PARTIAL DIFFERENTIAL EQUATIONS 101 6 4 10 The Fourth Order Runge Kutta Method in the Interaction Picture One way of reducing the memory overhead on the RK4 algorithm is to move into an inter action picture exactly as was performed
217. r this gives one the opportunity to choose whether or not to sample the initial point of the field and generate output moments such as mean and standard error before it is propagated To do this we put a sequence of 1 s or 0 s within the lt samples gt element We must put as many 1 s and 0 s as there are moment groups defined in the lt output gt tag to be discussed later In our case we will only be using one output moment group here and we do want to sample the initial field for the moments hence we add the code lt samples gt 1 lt samples gt to our script The simulation at this stage looks like cmi ee SO O lt simulation gt lt name gt lorenz lt name gt lt author gt Paul Cochrane lt author gt lt description gt Lorenz attractor example simulation Adapted from the example in Numerical methods for physics by Alejandro L Garcia page 78 1st ed lt description gt lt prop_dim gt t lt prop_dim gt lt globals gt lt CDATA const double sigma 10 0 const double b 8 0 3 0 const double r 28 0 const double xo 1 0 initial conditions const double yo 1 0 initial conditions 1 1 A BASIC SIMULATION 9 const double zo 20 0 initial conditions gt lt globals gt lt field gt lt name gt main lt name gt lt samples gt 1 lt samples gt lt yet more xmds code to come gt lt field gt lt simulation gt Now we have to initialise the field In othe
218. r words we have to define x t 0 y t 0 and z t 0 xmds makes this easy by allowing you to define the field as a vector written in terms of the dimensions of the field It is possible to define other vectors that are part of the field but the vector of the field that we are integrating is the main vector and this we name funnily enough main This naming is compulsory since it is possible to have more than one vector named within a field however there must be exactly one vector named main We also need to specify the data type of our vector the names of the components of the vector and we need to define how the vector should be calculated using C code in a CDATA section For the case we are considering here the xmds code would be lt vector gt lt name gt main lt name gt lt type gt double lt type gt lt components gt x y z lt components gt lt CDATA X XO y yo Z ZO gt lt vector gt So as we can see our vector is called main it is of type double and the names of its com ponents are x y and z The CDATA section gives the C code version of what Equation describes This code completes the lt field gt element and we are left with the following code listing lt xml version 1 0 gt lt simulation gt lt name gt lorenz lt name gt lt author gt Paul Cochrane lt author gt lt description gt Lorenz attractor example simulation Adapted from the example in Numerical method
219. rdinarily passed to the C compiler XMDS_LIBS the libraries and library paths necessary to build the simulation Again for those familiar with make this is the same as the LIBS variable XMDS_INCLUDES the include paths and files for the C C compiler to look for when compiling the simulation MPI options MPICC the C C compiler used by the local MPI implementation to compile the simulation Often this is just mpicc but on systems such as the APAC super computer in Canberra Australia this actually is just cc MPICCFLAGS the CFLAGS variable to be passed to the MPI C C compiler MPILIBS the extra libraries necessary to compile a simulation for use with MPI For instance on cluster systems an MPI implementation such as LAM will be used In this case the extra libraries necessary to compile a simulation will be something like lmpi llam FFTW options FFTW_LIBS the libraries and library paths specific to your fftw installation FFTW MPI LIBS the libraries and library paths specific to your fftw installation neces sary so that you can use fftw with MPI Warning only use fftw with a supercom puter To perform Fourier transforms in parallel a lot of communication between the nodes is necessary hence this is only worthwhile on a supercomputer with a high speed network connection between nodes FFTW3_LIBS the libraries and library paths specific to your installation of fftw version 3 FFTW3_THREADLIBS the additional
220. rences file here is an example xmds prefs xmds preferences file XMDS_CC gcc this here is a comment XMDS_CFLAGS pthread 03 ffast math funroll all loops XMDS_LIBS 1stdc 1m lxmds L home cochrane bin XMDS_INCLUDES I home cochrane bin THREADLIBS 1fftw_threads 130 FUNCTIONALITY The use of preferences can be switched on and off with the use of the lt use_prefs gt tag which should be located after the lt simulation gt element near where lt error_check gt and friends live This is a boolean option so to switch preferences off one should use no The default is yes but if you don t have any preferences or you haven t specified them all then the default values that were decided when xmds was built are used for anything not specified 8 2 Syntax summary The xml elements used can be divided into two main categories those that may contain code and or other elements and those that merely contain variable assignments The former can perhaps be called structural elements and the latter perhaps assignment elements These are summarised in Tables 8 2 and 8 3 respectively Element name Used in Reg May contain lt simulation gt top level yes lt name gt prop dim error check stochastic paths noises lt argv gt lt globals gt field sequence output lt binary_output gt use double u
221. ressing the issue of standards and reproducibility within computa tional physics The results of one group must be able to be independently verified by another group and should one group cease working in a particular area it should be not unnecessarily difficult for another group to pick up where they left off However code writing is inherently individualistic and does not naturally lend itself to these properties Further the range of problems that encompass every physical phenomenon known to our species is not small and if every problem had to be individually coded into a detailed computer program checked and debugged well the future would look bleak What is needed is to be able to specify the generalities of the problem without having to continually reprogram the detailed mechanisms of solving it This approach is known as high 79 80 INTRODUCTION level programming and there are a variety of both commercial and public license packages based on this approach Unfortunately to achieve the necessary flexibility these packages tend to be interpreters and thus they are seldom efficient at calculating the solution And for some problems those that take fast computers days and weeks to solve efficiency is paramount The only way to obtain efficiency is to directly compile a well written low level program letting the compiler do the optimisation which is what compilers are very good at What many mathematicians and scientists
222. rocessing required 146 WORKED EXAMPLE NLSE XMDS Consider firstly the lt sampling gt element Firstly a lt fourier_space gt assignment is com pulsory which contains as many yes or no entries as there are transverse dimensions so that the field may be transformed accordingly prior to sampling In this example the field s t dimension remains in normal space Next xmds looks for a lt lattice gt assignment containing the lattice on which to sample on for to sample the entire field may produce much more data than is necessary to obtain a good plot In this assignment there should be as many integers as there are transverse dimensions Each integer should be one of the following values or else xmds will exit with an appropriate error message e 0 A lattice integer of zero requests xmds to integrate the moments over this dimension This will cause the output field to no longer be a function of this transverse dimension e 1 A lattice integer of 1 requests xmds to sample the moments on a cross sectional slice of this dimension also causing the output field to loose this transverse dimension If this dimension is in normal space then xmds will extract the slice at the middle lattice point point number N 2 1 using integer division otherwise xmds will extract the slice at the zero momentum point k 0 e A factor of the main field s corresponding lattice value which is greater than one In this case if it is sam
223. rom A L Garcia Numerical Methods in Physics 1994 lt description gt lt prop_dim gt t lt prop_dim gt LE memo O mala prOopeyemlola Claim gt lt error_check gt yes lt error_check gt S COEUR LO Yes 2 lt use_wisdom gt yes lt use_wisdom gt lt I grama 10 mo gt lt benchmark gt yes lt benchmark gt CES dsrteaules o mo gt lt Global var ables fon the simulation gt gt gt lt globals gt lt CDATA const double v 1 0 const double x0 0 0 const double sigma 0 1 double k M_PI sigma gt lt globals gt Field to be integrated over gt lt field gt lt name gt main lt name gt lt dimensions gt x lt dimensions gt lt transverse dims gt lt lattice gt 50 lt lattice gt lt WMO Pus xo Caria chim gt lt domains gt 0 5 0 5 lt domains gt lt domain of each dim gt lt samples gt 1 lt samples gt ON sample ile DO E sidiam gt lt vector gt lt name gt main lt name gt lt type gt complex lt type gt lt l Gema WIS O VOCLOZ 32 lt components gt A lt components gt Si names Or components ee lt fourier_space gt no lt fourier_space gt lt def in k space gt lt CDATA A rcomplex cos k x x0 exp x x0 x x0 2 0 sigma sigma OOD ie 60 USING A TEMPLATE lt vector gt lt field gt SI hie sequence of integrations to perform gt lt sequence gt lt integrate gt algorith
224. rop_dim gt to t This time however we are performing a stochastic simulation so we need to keep the lt stochastic gt lt paths gt lt seed gt and lt noises gt tags Setting lt stochastic gt to yes the number of paths to 1024 leaving the lt seed gt setting as is and noting that we only have one noise term we get the following chunk of code describing the stochastic part of our simulation lt stochastic gt yes lt stochastic gt lt defaults to no gt lt these three tags only necessary when stochastic is yes gt lt paths gt 1024 lt paths gt IL no Of parias gt lt seed gt 1 2 lt seed gt lt I seeds fori AO eel gt lt noises gt 1 lt noises gt lt l WO gut MOLES gt We ll leave the use mpi tag where it is at present similarly with the rest of the tags down to the lt globals gt section We ve got three variables in this simulation omega the mean frequency gam the perturbation strength and zo the initial value of z The lt globals gt block comes out to this lt Global variables for the simulation gt lt globals gt lt CDATA double omega 0 mean frequency the word gamma is already used in many maths libraries hence we use gam here double gam 0 1 perturbation strength double zo 1 initial value of z 14 lt globals gt There are no transverse dimensions so we can get rid of the lt dimensions gt lt lattice gt and lt domai
225. s License If you cannot distribute so as to satisfy simultaneously your obligations under this License and any other pertinent obligations then as a consequence you may not distribute the Program at all For example if a patent license would not permit royalty free redistribution of the Program by all those who receive copies directly or indirectly through you then the only way you could satisfy both it and this License would be to refrain entirely from distribution of the Program If any portion of this section is held invalid or unenforceable under any particular circum stance the balance of the section is intended to apply and the section as a whole is intended to apply in other circumstances It is not the purpose of this section to induce you to infringe any patents or other property right claims or to contest validity of any such claims this section has the sole purpose of protecting the integrity of the free software distribution system which is implemented by public license practices Many people have made generous contributions to the wide range of software distributed through that system in reliance on consistent application of that system it is up to the author donor to decide if he or she is willing to distribute software through any other system and a licensee cannot impose that choice This section is intended to make thoroughly clear what is believed to be a consequence of the rest of this License 8 If the distribution and or us
226. s for physics by Alejandro L Garcia page 78 ist ed 10 STARTING FROM SCRATCH lt description gt lt prop_dim gt t lt prop_dim gt lt globals gt lt CDATA const double sigma 10 0 const double b 8 0 3 0 const double r 28 0 const double xo 1 0 initial conditions const double yo 1 0 initial conditions const double zo 20 0 initial conditions gt lt globals gt lt field gt lt name gt main lt name gt lt samples gt 1 lt samples gt lt vector gt lt name gt main lt name gt lt type gt double lt type gt lt components gt x y z lt components gt lt CDATA x xo M Re z zo Es lt vector gt lt field gt yet more xmds code to come gt lt simulation gt 1 1 4 The sequence element We have arrived at the stage where we can tell xmds how to actually perform the integration of the field To do this we use the lt sequence gt element The lt sequence gt element is usually used as a container for other elements specifically the lt integrate gt lt filter gt and other lt sequence gt elements The outermost lt sequence gt element is referred to as the parent lt sequence gt and the lt sequence gt s nested within that as the child lt sequence gt s This may sound a bit confusing but it is just a generalisation and a lot of the time you will be writing scripts with just the one lt sequence gt The lt sequence gt el
227. s gt lt CDATA temperature T gt lt sampling gt lt group gt lt output gt lt simulation gt Here is a link to the finished gzipped script file diffusion xmds gz on the xmds web site http www xmds org 1 2 7 Making the simulation and getting results Running through the sequence of events necessary to generate a simulation binary executable file running it and producing the results for Matlab Octave or Scilab you should see a 1 2 A MORE COMPLEX SIMULATION 29 sequence of events and output something like following Generating the simulation binary xmds diffusion xmds Running the simulation diffusion Making forward plan Making backward plan Beginning full step integration Sampled field for moment group 1 at t 0 000000e 00 Sampled field for moment group 1 at t 2 000000e 02 lt snip gt Sampled field for moment group 1 at t o O000 Oe SOn Sampled field for moment group 1 at t 1 000000e 00 maximum step error in moment group 1 was 9 909189e 09 The forward and backward plans are the fftw routines calculating the necessary Fourier transforms Generating the Matlab or Octave output xsil2graphics lorenz xsil Output file format defaulting to matlab Output file name defaulting to diffusion m Proccessing xsil data container 1 Writing data container 1 to file Generating the Scilab output xsil2graphics scilab lorenz xsil Output file name
228. s gt lt simulation gt 11 11 benchmark optional lt benchmark gt bool lt benchmark gt xmds 1 2 Contains boolean Subelement of Path to tag lt simulation gt lt benchmark gt Description Defaults to no Tells xmds whether or not to put timing code around the main section of code that part of the code excluding fftw creation and deletion as a way of benchmarking the simulation or at least giving an idea of how long the main section of code will take Example lt simulation gt lt benchmark gt yes lt benchmark gt lt simulation gt 11 12 binary output optional binary output bool lt binary_output gt xmds 1 2 OBSOLETE see output Contains boolean Subelement of Path to tag simulation binary output 11 13 USE_WISDOM 183 Description Defaults to no If set to yes xmds saves data file in binary format instead of the default ascii This tag is now obsolete Please specify the format attribute in the lt output gt tag For example lt output format binary gt Example lt simulation gt lt binary_output gt yes lt binary_output gt lt simulation gt 11 13 use wisdom optional lt use_wisdom gt bool lt use_wisdom gt xmds 1 2 Contains boolean Subelement of Path to tag simulation use wisdom Description Defaults to no If yes the simulation will use fftw s wisdom feature The simulation will look in the user s xmds wisdom directory
229. s listed in the lt dimensions gt lt components gt lt operator_names gt and lt moments gt assignments It also applies to user de fined variables where C code is allowed 8 3 C CODING WITHIN ELEMENTS 133 Function name Math Argument s Result abs n In integer integer fabs x x real real exp x e real real log x ln x real real sqrt x yT real real pow x y r real real real sin x sin x real real asin x sin x real real cos x cos x real real acos x cos l x real real tan x tan x real real atan x tan x real real cot x cotan x real real complex x y r real real complex rcomplex x y z real real complex pcomplex x y rey real real complex conj z or z z complex complex c_exp z ez complex complex c_log z In z complex complex c_sqrt z VE complex complex real z real z complex real imag z imag z complex real mod z Iz complex real arg z imag ln z complex real mod2 z z complex real Table 8 4 Commonly used functions Globals should be declared as asi J sj M const double reali 23 45 const complex compi complex 0 4 0 2 double r sqrt x x y y complex z pcomplex r M_PI 2 As can be seen there are three main variable types available integers reals and complex best declared as long double and complex respectively The range of mathematical fu
230. s when the integrate code is to be executed which will usually need to be after the moment groups are calculated 10 6 highdim vector version xmds SEI Ss O il 10 p lt simulation gt lt name gt highdim lt name gt Global system parameters and functionality gt lt prop_dim gt t lt prop_dim gt lt error_check gt yes lt error_check gt lt use_mpi gt yes lt use_mpi gt lt use_wisdom gt yes lt use_wisdom gt 168 MORE EXAMPLES lt benchmark gt yes lt benchmark gt lt Global variables for the simulation gt lt globals gt lt CDATA const double noise 0 0 const double hbar 1 05500000000e 34 const double M 1 409539200000000e 25 const double omegax 0 58976353090742 const double omegay 0 58976353090742 const double omegaz 0 58976353090742 30 const double U11 2 974797272874263e 51 const double U13 1 417820412490823e 50 const double U33 2 974797272874263e 51 const double inum 1 0e6 const double Uoh11 U11 hbar const double Uohi3 Ui13 hbar const double Uoh33 U33 hbar const double mu pow 15 inum U11 omegax omegay omegaz M PI 4 0 4 pow M 0 6 2 const double delta 1 0e9 const double F 2 0e 2 const double g sqrt Uohi1 2 0 delta const double lossii 1 0e 2 const double lossi2 1 6e 22 const double loss31 1 0e 2 const double loss32 1 6e 22 const double loss132 8 0e 17 const double chi F g delta const double biggamma g g delta 2 const double gam
231. s will search for a preferences file This file is called xmds prefs and can reside either in the user s HOME xmds directory or the directory local to the simulation script being processed xmds searches for the file in the user s home directory first and then in the current directory 2 8 2 Setting preferences Preferences are set in the xmds prefs file by using key value pairs delimited by an equals sign Therefore in general the format is key value Spaces around the equals sign are ignored such that the following are all equivalent key value key value key value key value The hash character is used for comments anything after and including the are ignored when parsing the options So one can make what is happening much clearer what the flags are to do For example 48 EXTRA AND ADVANCED FEATURES this is my funky preferences file option some_funky_value some_other_funky_value other_option also_quite_funky_variable isn t this variable funky the first line will be ignored and the text including and after the on the third line will be ignored 2 8 3 What are the options The options for building simulations that xmds accepts are e xmds options XMDS_CC the C C compiler xmds will use Typical options include cc gcc g XMDS_CFLAGS the flags passed to the C C compiler For those who are using the make utility then these flags are the same as the CFLAGS variable o
232. sampling gt lt fourier_space gt lt fourier_space gt lt sample in k space gt lt lattice gt lt lattice gt lt l mo points to sample gt lt moments gt lt moments gt lt names of moments gt lt CDATA gt lt sampling gt lt group gt lt output gt lt simulation gt 3 3 Ripping it to bits Now all we need to do is to work out what we need and don t need and to throw the relevant numbers and equations into the relevant places If the pace is too high you might like to go back and re read Chapter 1 to make sure you got it all 3 3 1 Global system parameters and functionality Looking at the block of code that sits between the start lt simulation gt tag and the lt globals gt element we can add the following things the name of the simulation is advection the author is Paul Cochrane and the description can be something like lt description gt Solves the one dimensional advection equation for an initial cosine modulated Gaussian pulse Adapted from A L Garcia Numerical Methods in Physics 1994 lt description gt The propagation dimension is in time so we set lt prop_dim gt to t the simulation isn t stochastic so we can leave the code where it is or remove it since it s set to the default setting anyway We ll remove it in this case as we don t gain anything and it s extra code floating around that we don t need Since this isn t a stochastic simulation t
233. sdom where the lt rank gt is the MPI process rank number and stops name conflicts when doing parallel simulations There are two places that xmds can store wisdom files in the user s xmds wisdom directory or in the directory local to the simulation The former is used if the xmds wisdom directory exists and the latter is used if not Running the atomlaser simulation with wisdom turned on we get the following output for the first run Performing fftw calculations Making forward plan Making backward plan Keeping accumulated wisdom Finished fftw calculations Beginning full step integration Sampled field for moment group 1 at t 0 000000e 00 Sampled field for moment group 1 at t 2 500000e 09 Sampled field for moment group 1 at t 5 000000e 09 Sampled field for moment group 1 at t 7 500000e 09 Sampled field for moment group 1 at t 1 000000e 08 Beginning half step integration Sampled field for moment group 1 at t 0 000000e 00 Sampled field for moment group 1 at t 2 500000e 09 Sampled field for moment group 1 at t 5 000000e 09 Sampled field for moment group 1 at t 7 500000e 09 Sampled field for moment group 1 at t 1 000000e 08 maximum step error in moment group 1 was 3 802825e 11 Time elapsed for simulation is 11 seconds and then for the second run Performing fftw calculations Standing upon the shoulders of giants Importing wisdom Making forward plan Making
234. se wisdom lt benchmark gt use prefs gt lt argv gt threads fftw version lt globals gt simulation no C code lt argv gt lt simulation gt no lt arg gt lt arg gt lt argv gt yes lt name gt lt type gt lt default_value gt lt field gt lt simulation gt yes lt name gt lt dimensions gt lt lattice gt lt domains gt lt samples gt lt vector gt lt vector gt lt field gt yes lt name gt lt type gt lt components gt lt vectors gt lt filename gt lt fourier_space gt and C code lt sequence gt lt simulation gt yes lt integrate gt lt filter gt lt sequence gt lt sequence gt no lt cycles gt lt integrate gt lt filter gt lt sequence gt lt integrate gt lt sequence gt no lt algorithm gt lt interval gt lt lattice gt lt samples gt lt k_operators gt lt vectors gt lt cross_propagation gt lt iterations gt lt moment_group gt lt tolerance gt lt max_iterations gt lt min_time_step gt lt cutoff gt lt smallmemory gt lt no_noise gt lt halt_non_finite gt and C code lt k_operators gt lt integrate gt no lt constant gt lt operator_names gt and C Code lt cross_propagation gt lt integrate gt no lt prop_dim gt vectors and C Code lt filter gt lt sequence gt no lt fourier_space gt lt vectors gt lt functions gt lt moment_group gt and C Code lt output gt lt simulation gt yes
235. seful if one wants to vary pa rameters within the simulation between simulation runs Hence one can write a Perl Python shell script to run the simulation over the various parameters with out having to rewrite an xmds script recompile and then rerun Example lt simulation gt lt argv gt ME mdi lt argv gt lt simulation gt 11 20 1 arg required lt arg gt xmds tags lt arg gt xmds 1 2 Contains lt default value gt Subelement of Path to tag lt simulation gt lt argv gt lt arg gt Description Container for the tags describing a given command line argument 188 LANGUAGE REFERENCE Example lt simulation gt lt argv gt lt arg gt SU sas vela gt gt lt arg gt lt argv gt lt simulation gt 11 20 2 name arg required lt name gt char lt name gt xmds 1 2 Contains string Subelement of Path to tag lt simulation gt lt argv gt lt arg gt lt name gt Description The name of the command line argument This tag is used to create a flag to specify the value of the variable entered at the command line Two forms are accepted a short form and a long form each conforming to the GNU getopt conventions For instance for a variable named nconst the flag used at the command line will have a short form of n and a long form of nconst However if another variable has been chosen beginning with the letter n then n is no longer unique and to make the
236. sequence gt lt integrate gt lt algorithm gt RK4EX lt algorithm gt lt interval gt 1 lt interval gt lt lattice gt 1000 lt lattice gt lt samples gt 50 lt samples gt lt k_operators gt lt constant gt yes lt constant gt lt operator_names gt L lt operator_names gt lt CDATA L kappa kx kx 1 gt lt k_operators gt 26 STARTING FROM SCRATCH lt CDATA drid de Wes lt integrate gt lt sequence gt 1 2 5 The output element The lt output gt element is much like that discussed in Section We need to specify a lt filename gt tag just so that we are documenting everything nicely and then lt group gt and lt sampling gt elements to describe the rest of the information xmds needs to properly sample the data and save it out to file There are two new tags that are needed because we have a transverse dimension to worry about these are lt fourier_space gt and lt lattice gt The lt fourier_space gt tag is necessary to tell xmds if the sampling of the output is to be performed in Fourier space xmds expects to see a space separated list of yes or no values for each of the transverse dimensions and for the situation here we don t want the output sampled in Fourier space and we only have one transverse dimension so we set the lt fourier_space gt tag to no The lt lattice gt tag tells xmds how finely the transverse dimensions should be sampled and in general expects to see a space separated lis
237. simulation gt lt sequence gt lt integrate gt lt k_operators gt lt operator_names gt L lt operator_names gt lt k_operators gt lt integrate gt lt sequence gt lt simulation gt 11 22 2 13 moment group integrate optional moment group xmds tags lt moment group Contains lt integrate_dimension gt CDATA Subelement of Path to tag lt simulation gt lt sequence gt lt integrate gt lt moment_group gt Description Defines and calculates a number or vector integrated through zero or more of the transverse dimensions of the problem 208 LANGUAGE REFERENCE Example lt simulation gt lt sequence gt lt integrate gt lt moment_group gt lt moments gt joe is great lt moments gt lt integrate_dimension gt no yes lt integrate_dimension gt lt CDATAT joe psi psi is phi phi great theta gt lt moment_group gt lt integrate gt lt sequence gt lt simulation gt 11 22 2 14 functions integrate optional lt functions gt CDATA lt functions gt Contains CDATA Subelement of Path to tag lt simulation gt lt sequence gt lt integrate gt lt functions gt Description This is the best place to define any functions that do not depend on the transverse dimensions Example lt simulation gt lt sequence gt lt integrate gt lt functions gt lt CDATA Some C code I gt lt functions gt lt integr
238. split lines of code over multiple lines if necessary 24 STARTING FROM SCRATCH and the M_PI variable The rcomplex function is one of a set of utility functions added as part of xmds to allow users to define complex variables The syntax of rcomplex is rcomplex x y where x and y are real variables representing the real and imaginary parts of the complex number respectively This explains one line of the CDATA block reads merely 0 0 The lonely comma is just separating the two arguments to rcomplex the 0 0 is the imaginary part of the variable we are defining which is real but of complex type and the closing parenthesis and semicolon just finish off the function call syntax Notice that the variable T is assigned to a quantity which is defined over multiple lines Although this may seem strange to some for those familiar with C language rules will know that all the C compiler is looking for is the semicolon as the character denoting the end of the expression Therefore it is possible to split equations over many lines to break up complicated expressions or to highlight certain parts of the expression that may be important The M PI variable is an automatically set variable by xmds that is the value of 7 i e 3 14159 Note also that one can use any of the standard mathematical functions defined in C C such as sqrt log etc This completes the discussion of the lt field gt element which is now I Field to be in
239. sses the input file and populates a node tree based on the Document Object Model IT This parser was written by the author but there is little need to describe it in detail here This node tree is then passed to a simulation class object which extracts its own relevant data and in turn creates child objects to process the relevant sub trees The class structure of xmds is shown in Figure Once the input file has been processed a dry run of the main sequence tree is called in order to evaluate the total number of samples requested for each output moment group Finally the output code is generated This is accomplished with a tree walking tech nique To begin with the simulation object writes any include statements that are necessary Then it writes any define statements particular to itself and then calls each of its child objects to do the same This process is continued down the tree until every object in the tree has written its define statements Starting again at the simulation element the process is repeated to write the global variables the routine prototypes and finally the routines themselves This tree walking is performed by the Element class with the vir tual functions aptly titled writeDefines writeGlobals writePrototypes and writeRoutines These routines are then overridden in the derived classes to write the code specific to each class As mentioned in the previous section one of the primary areas of difficu
240. st be binary not ascii Also the moment group number of the XSIL file that will be used for initialisation must be specified with the moment_group attribute of the lt filename gt tag if there is more than one moment group in the XSIL file if there is only one moment group in the XSIL file as is the case for XSIL files generated from a lt breakpoint gt tag then this attribute can be omitted If the vector is of type double then the variables of the vector are initialised from the output moment group variables of the same name but suffixed with an R If the vector is of type complex then the real and imaginary components of each variable are initialised by the values of the output moment group variables of the same name but with a suffix of R for the real component and T for the imaginary component For example the complex variables x and y would be initialised by the output moment group variables xR xI yR and yI and you would use the following code in your output tag to create these variables lt output gt lt group gt lt sampling gt lt moments gt xR xl yR yl lt moments gt lt CDATA SR HORES x x im 40 EXTRA AND ADVANCED FEATURES 13 y mer Wil Yo imi ie lt sampling gt lt group gt lt output gt Not every variable in a vector need be present in the moment group of the XSIL file as any variable that is not present is automatically initialised to zero or by a CDATA secti
241. st to see how long the main body of code takes you can use the lt benchmark gt tag This tag is a boolean which by default is no but when set to yes it tells xmds to insert timing code around the main code block excluding the fftw plan creation and deletion steps this is because these steps are not in general indicative of how long the simulation will run especially when scaled up to long simulation times To use this feature just put the code lt benchmark gt yes lt benchmark gt into your simulation script in the global simulation and functionality section namely before the lt globals gt tag The simulation will then report at the end of its run how long it took An example of this is again using the atomlaser simulation Making forward plan Making backward plan Beginning full step integration Sampled field for moment group 1 at t 0 000000e 00 Sampled field for moment group 1 at t 2 500000e 09 Sampled field for moment group 1 at t 5 000000e 09 Sampled field for moment group 1 at t 7 500000e 09 Sampled field for moment group 1 at t 1 000000e 08 Beginning half step integration Sampled field for moment group 1 at t 0 000000e 00 Sampled field for moment group 1 at t 2 500000e 09 Sampled field for moment group 1 at t 5 000000e 09 Sampled field for moment group 1 at t 7 500000e 09 Sampled field for moment group 1 at t 1 000000e 08 maximum step error in moment group 1
242. t CDATA dffi_dz Lap1 ff1 i ff2 sh damping ff1 dff2 dz Lapi ff2 i ffi sh damping ff2 dsh_dz Lap2 sh i ff1 ff2 damping sh gt lt integrate gt lt sequence gt lt output gt lt group gt lt sampling gt lt fourier_space gt no no lt fourier_space gt lt lattice gt 50 O lt lattice gt lt moments gt pow_dens lt moments gt lt CDATA pow dens ffl ffl ff2 ff2 2 sh sh Jie lt sampling gt lt group gt lt group gt lt sampling gt lt fourier_space gt no no lt fourier_space gt lt lattice gt 0 0 lt lattice gt lt moments gt etot lt moments gt lt CDATA etot ffi ffi ff2 ff2 2 sh sh Jal lt sampling gt lt group gt lt output gt lt simulation gt This simulation describes how to solve for the evolution of three bosonic fields governed by a non degenerate parametric interaction as described in Equations and 10 2 Perhaps not surprisingly the majority of the above script is to generate the particular ini tialisation conditions and to define a number of output moments the portion concerned with actually implementing these equations is not that great Oo 9 9 x tato 1 6 ek s j 12 10 1 o o ci IUE Ga oda 10 2 152 MORE EXAMPLES The main difference between this simulation and the nlse xmds simulation is that it now has a three component field which has two transverse dimensions Also two
243. t lt field gt In the lt sequence gt block we set the algorithm to SIEX because a semi implicit algorithm does a better job of integrating stochastic equations than the fourth order Runge Kutta We ll explicitly set the number of iterations to the default of 3 we might want to change it in the future we ll integrate for 10 seconds so set interval to 10 We ll use 1000 for the lattice size and take 100 samples of it in the output moment group There aren t any k operators here so we can delete that entire section and the lt vectors gt tag The interesting bit is where we write in the differential equation to be solved which is really easy to do and is written in the CDATA block as lt CDATA dz_dt i omega z i sqrt 2 0 gam n 1 z 11 gt where we ve expanded out the brackets of the Kubo oscillator equation Equation 4 1 Notice that the value of z mentioned in Equation has been replaced by the variable n 1 This is because when we tell xmds that we have noises it automatically creates some noise variables for us In the present case we only have one noise so there is only one variable namely n 1 However in the general case where we could have say m noises we would have the noise variables n 1 n 2 n m Putting these pieces together gives a sequence block that looks like this gt gt lhe sequence ot integrations to perform gt lt sequence gt lt integrate gt lt algorithm gt SIEX lt algorithm
244. t lt lattice gt 4 8 16 lt lattice gt lt moments gt atomsr moleculesr atomsi moleculesi lt moments gt lt CDATA atomsr phi 1 moleculesr phi 3 atomsi ixphi 5 moleculesi ixphi 6 Tp lt sampling gt lt group gt lt output gt lt simulation gt This simulation is identical in function to the highdim xmds example above but describes the fields as an array of components rather than a list This notation may be very valuable 10 6 HIGHDIM_VECTOR_VERSION XMDS 173 when the numbers of component get very large and the equations can be easily described in terms of the index WARNING There is no bounds checking on the index of your field so be careful when writing your equations in this form When using this notation if XMDS needs to calculate the k space operator of any of the components of an array all of them are calculated This makes for example high dim_vector_version xmds slower than the old version IMPORTANT For interaction picture algorithms if a k space operator is applied to any component of a vector then it is applied to ALL OF THEM This means that high dim_vector_version xmds only solves the correct equations when used with an EX algorithm 174 MORE EXAMPLES Part IV Reference Manual 175 Language Reference 11 1 simulation required lt simulation gt xmds tags lt simulation gt Contains prop dim lt gTobaTs gt usempi IPI Method lt field Ksequence K
245. t and introduce some new xmds tags but still with the theme that we are writing this all from a clean slate Hopefully you will be able to see the other more powerful abilities of xmds and be able to start writing your own simulations 1 2 A more complex simulation In this section I ll introduce how to use Fourier space in your simulations and the extra tags required and explain why it is sometimes easier to perform part of the calculation in Fourier space and then transform back to position space To illustrate these extensions to what we already know from Section we ll look at solving the one dimensional diffusion equation 2 Oa z t _ a x t 14 Ot Ox where a z t is the field to be evolved by the differential equation and is a function of time t and space z and amp is a constant describing how quickly the solution diffuses An example of an application of the diffusion equation is for modelling the diffusion of some initial distribution of temperature in a metal rod over time the temperature distribution will flow from areas of higher temperature to areas of lower temperature eventually achieving a uniform distribution over the entire rod 1 2 A MORE COMPLEX SIMULATION 21 So why do we use Fourier space when solving this differential equation The main reason is that it s a lot easier to calculate some of the differentials in Fourier space than it is in position space It turns out that if one transforms position sp
246. t fourier_space gt lt CDATA const double wO hwhm sqrt 2 log 2 const double amp sqrt energy w0 sqrt M_PI 2 phi pcomplex amprexp t t w0 w0 velx t Wile lt vector gt The lt name gt assignment is compulsory here as it is quite common to employ more than one vector The field must have at least one vector called main as it is this vector that xmds forward evolves in the main propagation dimension The lt type gt assignment is optional It will default to complex unless specifically stated as double which is the standard data type in C for floating point variables Note that vectors of type double cannot be Fourier transformed without some arbitrary definition of how this is done Thus xmds requires that double vectors always remain in the Fourier 9 5 THE SEQUENCE ELEMENT 141 space in which they were initialised A request to access them in some other Fourier space will result in an error The lt components gt assignment is compulsory It identifies by name the components of the vector Here there is only one component but more may be defined simply by separating them with spaces for example lt components gt phi theta lt components gt The component names must be valid as far as the C language is concerned see Section 8 3 Finally the initialisation of the vector must be defined This can be done in two ways Either by means of a lt filename gt assignment or by sup
247. t all of the data is contained within the XSIL file However with binary output the data is pointed to by the ascii XSIL file where the actual data is distributed among different files labelled according to simulation name and output moment group In general the filenames are combined like this simulation name mg moment group number dat For instance if the simulation name is nlse then the output binary file for the first moment group will be called nlsemg1 dat In order that old output data remains meaningful xmds copies the input simulation script to the name of the output file overwriting it if it already exists and then the output data is inserted as lt XSIL gt elements before the closing lt simulation gt tag This way the input parameters that generated the output data remain known Next xmds looks for the lt group gt elements contained within of which there must be one or more The lt group gt element contains moments that are sampled on a common output lattice in a common Fourier space If necessary the group may then be post processed after all segments have been performed for example to calculate moments with the propagation dimension itself in Fourier space Thus xmds looks for a compulsory lt sampling gt element within the lt group gt element that defines how to sample for the moments It then looks for an optional lt post_propagation gt element within the same lt group gt element that defines any post p
248. t area of the field to have no effect on the result and when periodic boundary conditions are used these are in fact enforced with xmds damping must be applied near the boundaries so as to absorb any such radiation The lattice used in each transverse dimension must be sufficiently fine to be able to exhibit the details yet not so fine as to cause memory problems or cause the solution to take forever to compute Further still field details that are significantly smaller than those present at initialisation may evolve during the simulation and the evolution of such detail will ultimately be affected by the finesse of the transverse lattice Therefore simply reducing the time step alone is not a complete guarantee of accuracy The results must be inspected to ensure that the transverse lattice was fine enough to contain the detail If unsure then re run the simulation with a finer lattice as the evolution may be attempting to generate singularities Usually meaningful results only come from stochastic integration by averaging moments of the fields over large numbers of paths This leads to a second form of error for stochastic equations due to the standard error in the mean This is known as the sampling error and scales as the inverse square root of the number of paths taken 114 NUMERICAL MODELLING THEORY Part III User Manual 115 Development and Program Structure It only takes a few items of information to define a numerical simu
249. t of values telling it how many samples in the particular transverse dimension to take Here we just set lt lattice gt to 50 The moments we are interested in sampling is of the temperature so we specify a lt moments gt tag value of temperature and give the code to calculate this moment in a CDATA block as lt CDATA temperature T gt All of this information gives the lt output gt element to be S The output to generate gt lt output gt lt filename gt diffusion xsil lt filename gt lt group gt lt sampling gt lt fourier_space gt no lt fourier_space gt lt lattice gt 50 lt lattice gt lt moments gt temperature lt moments gt lt CDATA temperature T gt lt sampling gt lt group gt lt output gt 1 2 A MORE COMPLEX SIMULATION 27 1 2 6 The final script We now have a complete script And this is in its entirety lt ame rS ion i m S locale SLELO Dakara MSO ACTE O gt lt simulation gt Global system parameters and functionality gt lt name gt diffusion lt name gt lt author gt Paul Cochrane lt author gt lt description gt Solves the one dimensional diffusion equation for an initial Gaussian pulse Adapted from A L Garcia Numerical Methods in Physics 1994 lt description gt lt prop_dim gt t lt prop_dim gt lt l Global variables fon the simulation gt gt lt globals gt lt CDATA const double kappa 0 1 diffusion coe
250. t s Contains integer 198 LANGUAGE REFERENCE Subelement of Path to tag lt simulation gt lt sequence gt lt cycles gt Description The number of times to perform a given sequence as a subsequence of the main sequence Defaults to 1 Example lt simulation gt lt sequence gt lt cycles gt 3 lt cycles gt lt sequence gt lt simulation gt 11 22 2 integrate optional lt integrate gt xmds tags lt integrate gt Contains FREE lt smallmemory gt lt cutoff gt lt halt_non_finite gt lt lattice gt lt samples gt lt k_operators gt lt noment_group gt functions gt lt vectors gt CDATA Subelement of Path to tag lt simulation gt lt sequence gt lt integrate gt Description Container element holding the tags that describe how the integration should take place Example lt simulation gt lt sequence gt lt integrate gt Slc made ages gt lt integrate gt lt sequence gt lt simulation gt 11 22 2 1 algorithm optional lt algorithm gt string algorithmName lt algorithm gt Contains string one of RK4EX RKAIP ARK45EX ARK45IP SIEX or SIIP Subelement of 11 22 SEQUENCE 199 Path to tag lt simulation gt lt sequence gt lt integrate gt lt algorithm gt Description The algorithm to use when integrating the equations Defaults to SIEX if lt stochastic gt is yes or RK4EX if lt stochastic gt is no The si
251. t up scilab and at its command prompt run the command gt exec lorenz sci gt temp_di zeros 1 201 gt t_1 zeros 1 201 x ut 1 zeros 1 201 yOut 1 zeros 1 201 z ut 1 zeros 1 201 error x ut 1 zeros 1 201 error y ut 1 zeros 1 201 error z ut 1 zeros 1 201 lorenzi fscanfMat lorenzi dat Error Info buffer is too small too many columns in your file 7 gt temp_d1 lorenz1 1 x ut 1 lorenz1 2 gt yOut_1 lorenz1 3 z 0ut 1 lorenz1 4 gt error_xOut_1i lorenz1 5 error y ut 1 lorenzi 6 gt error_zOut_1i lorenz1 7 gt t_1 temp di 2 clear lorenzi temp_di to load the data into Scilab you can safely ignore the warning and then to obtain a graphical output of the data run the command 20 STARTING FROM SCRATCH param3d xOut 1 y0ut_1 z0ut_1 which should give something along the lines of that in Figure 16 0 FIGURE 1 2 Three dimensional plot in Scilab of the trajectories of a Lorenz attractor Pa rameters used were c 10 b 8 3 r 28 with initial conditions of zo 1 0 yo 1 0 and 20 20 0 As we can see from both of these figures that we get the usual strange attractor butterfly shape Now that we ve spent a lot of time going over the very basics of writing a simulation from scratch we now speed up a bi
252. t uses the values of the field from the other lattice points There are two main methods for doing this 88 NUMERICAL MODELLING THEORY Firstly there is the matrix method If we let a represent the discretised values of a x at the lattice points x then we may evaluate the first order derivative like so EN ETETE E Ox Zk 1 Zk 1 M ra k F 1 n 6 29 Note the symmetry of the weightings about the point where the derivative is being evaluated causing the method to be space centered Note also that unless periodic boundary conditions are specified it is not not possible to implement a space centered method at the ends of the lattice The best that can be done is to implement a non space centered method at the ends and insist that the derivatives vanish at the lattice boundaries for the technique to be valid M may be expressed as a tri diagonal square matrix which for periodic boundary conditions is 0 1 0 0 e ce 1 1 0 1 O 1 0 0 1 0 we 6 30 z s amp d d 0 0 1 0 1 1 6 4 D oL 0 And similarly for evaluating the second derivative we might use 2 1 0 0 1 1 1 0 0 1 0 1 2 1 0 M ne arn mee 6 31 2 Colo is a ed u 0 0 1 1 lo 0 0 1 2 These methods evaluate the derivatives using only the adjacent points in space but in principle it is possible to use progressively more distant neighbours to produce an estimate of each derivative which takes i
253. t vector gt lt filename format xsil moment_group N geometry matching mode strict loose gt Si enter thel filename here gt lt filename gt lt fourier_space gt lt yes no gt lt fourier_space gt lt CDATA optional CDATA code Jie lt vector gt lt field gt where strict mode is the default geometry matching mode 2 6 2 Input data layout for ASCII and binary formats We now know the syntax of how to tell xmds that we want to input data from file we just now need to organise the data that we are going to input into the layout that xmds expects to see it Let s see how this works by considering a simple example Imagine we have three input vectors that we want to initialise with double precision data x y and z Their values are x 2 0 1 0 0 0 1 0 2 0 42 EXTRA AND ADVANCED FEATURES 5e 2 1e 3 1e 5 2e 4 7e 2 10 20 30 50 1e3 y z We can see that they are all 5 elements long this will equal the lt 1attice gt assignment and that they can contain numbers formatted in exponential notation We ll save this data into a file called input dat xmds expects this data to be ordered in a particular way which is related to the way the data is stored internally This order is an interlacing of the elements of each vector such that the first element of the first vector in this case x is expected as the first entry in the input file then the first element of the second vecto
254. te a 10 4 TLA XMDS 159 lt CDATA const double g 1 const double t0 1 gt lt globals gt lt field gt lt dimensions gt t lt dimensions gt lt lattice gt 100 lt lattice gt lt domains gt 10 15 lt domains gt lt samples gt 1 0 lt samples gt lt vector gt lt name gt main lt name gt lt type gt double lt type gt lt components gt E lt components gt lt CDATA E 2 t0 cosh t t0 gt lt vector gt lt vector gt lt name gt cross lt name gt lt type gt double lt type gt lt components gt P N lt components gt lt CDATA DEL 0 WN i gt lt vector gt lt field gt lt sequence gt lt integrate gt lt algorithm gt RK4EX lt algorithm gt lt interval gt 4 lt interval gt lt lattice gt 50 lt lattice gt lt samples gt 50 50 lt samples gt lt vectors gt main cross lt vectors gt lt CDATA dE_dz g P ike lt cross_propagation gt lt vectors gt cross lt vectors gt lt prop_dim gt t lt prop_dim gt lt CDATA dP_dt E N 160 MORE EXAMPLES dN_dt ExP I lt cross_propagation gt lt integrate gt lt sequence gt output group sampling lattice 50 lt lattice gt moments pow dens lt moments gt lt CDATA pow_dens ExE 1 1 gt lt sampling gt lt group gt lt group gt lt sampling gt lt vectors gt main cross lt vectors gt lt lattice
255. tegrate step a filter step or else a sub loop of further segments The key to xmds as explained in the introduction is that it is a code generator xmds requires the user to write their particular PDEs or ODEs as a few lines of C code which are 8 1 INSTALLING AND RUNNING xmds 125 Explicit picture Interaction picture Semi Implicit SIEX SIIP 4th Order Runge Kutta RK4EX RK4IP 4th 5th Order adaptive Runge Kutta ARK45EX ARK45IP 9th Order Runge Kutta RK9EX RK9IP 8th 9th Order adaptive Runge Kutta ARK89EX ARKS9IP Table 8 1 The Algorithm Matrix transplanted from the input script to the relevant points in the output code This technique proves itself well for efficient code that is extremely flexible A potential user should not be daunted at the prospect of having to learn basic C syntax the flexibility gained by this approach well justifies the effort In terms of output results it is usually desired to sample the field at various points throughout the sequence of operations and thus generate a sequence of samples Further when sampling often the raw complex value of the field at the point in question is not relevant it is a more general moment of field components and dimensions that is desired The examples provided in Chapters 9 and 10 will help explain Also the evolution of the field is usually solved on a lattice much finer than is necessary to make a good plot of the output Thus to sa
256. tegrated gt lt field gt lt name gt main lt name gt lt dimensions gt x lt dimensions gt lt lattice gt 100 lt lattice gt lt domains gt 1 1 lt domains gt lt samples gt 1 lt samples gt lt vector gt lt name gt main lt name gt lt type gt complex lt type gt lt components gt T lt components gt lt fourier_space gt no lt fourier_space gt lt CDATA T rcomplex exp x x0 x x0 2 0 sigma sigma sigma sqrt 2 0 M_PI 0 0 1 gt lt vector gt lt field gt 1 2 4 The sequence and integrate elements Using what we know from Section 1 1 we now use lt sequence gt and lt integrate gt elements to describe the guts of what xmds has to do We ll use here the RK4EX algorithm an interval 1 2 A MORE COMPLEX SIMULATION 25 of length 1 a lattice of 1000 and take 50 samples along the propagation direction Now because we re evolving the solution partially in Fourier space we need to define the operators that are going to be performing the evolution This is done with the lt k_operators gt element The reason why we call this the lt k_operators gt element is because Fourier space is often referred to as k space and position space as x space hence these operators are operating in k space and so they are k operators We next tell xmds that the k operators there is actually only one here are constant over the course of the simulation by setting the lt constant gt
257. tep defined in the simulation script via the lt lattice gt and lt interval gt assignments and then again at half of the full step size The maximum difference between the field values in each moment group is reported and gives an indication of the discretisation error in the simulation If one of the adaptive algorithms is chosen error checking means that the simulation is run a second time with one 16th of the specified tolerance For instance setting lt error_check gt to yes in the atomlaser simulation this is an example code in the examples directory of the xmds distribution we get the following output Making forward plan Making backward plan Beginning full step integration Sampled field for moment group 1 at t 0 000000e 00 Sampled field for moment group 1 at t 2 500000e 09 Sampled field for moment group 1 at t 5 000000e 09 Sampled field for moment group 1 at t 7 500000e 09 Sampled field for moment group 1 at t 1 000000e 08 Beginning half step integration Sampled field for moment group 1 at t 0 000000e 00 Sampled field for moment group 1 at t 2 500000e 09 Sampled field for moment group 1 at t 5 000000e 09 Sampled field for moment group 1 at t 7 500000e 09 33 34 EXTRA AND ADVANCED FEATURES Sampled field for moment group 1 at t 1 000000e 08 maximum step error in moment group 1 was 4 408207e 11 The error reported here is of the order of 107 and therefore we ca
258. ter code GCC 4 seems to have no problem There seems to be no problem using version 10 1 of Intel s compiler icc Note that icc sacrifices floating point compliance for speed by default but this can be changed using the fp model flag Remember you can change compilation flags by editing the preference file see section 2 8 for more information on preferences Of course you can also edit the compilation command printed by xmds and run it by hand Overall however there are no guarantees that halt_non finite will work when sac rificing accuracy for speed in floating point arithmetic Defaults to no Example lt simulation gt lt sequence gt lt integrate gt lt halt_non_finite gt yes lt halt_non_finite gt lt integrate gt lt sequence gt lt simulation gt 11 22 2 10 lattice integrate required lt lattice gt int lt lattice gt Contains integer Subelement of Path to tag lt simulation gt lt sequence gt lt integrate gt lt lattice gt Description The number of points to use over the integration interval i e over the number entered in the lt interval gt assignment Example lt simulation gt lt sequence gt lt integrate gt lt lattice gt 1000 lt lattice gt lt integrate gt 11 22 SEQUENCE 205 lt sequence gt lt simulation gt 11 22 2 11 samples integrate required samples int int lt samples gt Contains array of integers Sube
259. ter we will outline the tags necessary for xmds to actually parse a document and will implement a simple simulation involving these tags More complete documentation of the tags and what they do can be found in Chapter 11 The tags necessary for performing more advanced operations and for solving more complex problems will be introduced in later tutorial chapters e g see Chapter 4 1 1 A basic simulation 1 1 1 The simple beginnings of a simulation xmds is coded using XML the extensible markup language and as such each xmds script must be a properly formed XML document One of the main stipulations for an XML document to be properly formed is that it have the following line at the top of each document and so each xmds script must have this as its first line lt xml version 1 0 gt There are other stipulations but they don t really concern us too much at the moment If you re interested you can check out the World Wide Web Consortium W3C web site http www w3c org for more information 4 STARTING FROM SCRATCH Each xmds simulation is enclosed within a set of lt simulation gt tags Therefore for each simulation that you write from scratch the first few lines of code are going to be lt xml version 1 0 gt lt simulation gt AI yet more xnds code to come gt gt gt lt simulation gt It might be a good idea at this stage to save this code to file So for the purposes of the tutorial we
260. ternating steps hence this method is also known as a Split Step method or a Split Operator method The exact implementation depends on the main integration method employed as is detailed in the following Sections and As a general rule the interaction picture algorithms are faster and more stable than their Explicit picture counterparts but they are also more restrictive with regard to the allowable linear operators 6 4 9 The Semi Implicit method in the Interaction Picture The original derivation of this method was performed by Drummond 6 and we repeat it here since it is the most efficient and stable method for solving stochastic PDEs Consider now expressing the field a as the transform ai x x1 2 e oe b 29 x1 6 38 The time derivative of a then becomes aout 25 x2 2 k ai 7976 xo ie 6 39 100 NUMERICAL MODELLING THEORY since the operators are linear Equating this result with Equation 6 37 we get bx x e 72 2931 N x a 6 40 Now solving for b using the semi implicit method yields TA NEL i 0 b z y Xl x x4 6 41 h e 29 2 C 2 1 i G dle A xi a J ju 6 42 Finally if we define z z9 i then a and b become identical at the middle of the time step and we obtain p 4 3 b ax TA Gob 2 3 6 43 Hence the a and b fields are transformed from one to the other at the beginning and end of the time step using the linea
261. the CDATA block must conform to C C syntax rules otherwise the simulation won t be able to be compiled This is because the code within the CDATA block is inserted directly into the code for the output simulation Secondly the constants that we are specifying are declared to be double precision to xmds this is via the double keyword in the code above since they are continuous variables in the problem being solved If for instance we had some discrete quantity such as the number of particles in a system then we would specify the variable as being integer and would therefore use the int keyword to declare the variable as such Lastly the odd looking tags for the CDATA block must be written correctly i e opened with lt CDATA and closed with gt otherwise the file won t parse and xmds will give an error Our script is now lt xml version 1 0 gt lt simulation gt lt name gt lorenz lt name gt lt author gt Paul Cochrane lt author gt lt description gt Lorenz attractor example simulation Adapted from the example in Numerical methods for physics by Alejandro L Garcia page 78 1st ed lt description gt lt prop_dim gt t lt prop_dim gt lt globals gt lt CDATA const double sigma 10 0 const double b 8 0 3 0 const double r 28 0 const double xo 1 0 initial conditions const double yo 1 0 initial conditions const double zo 20 0 initial conditions 11 lt globals gt
262. tion Re capping the algorithms fell into two main categories 1 Explicit picture As detailed in Section 6 4 2 the transverse derivatives of the field are calculated using the Fourier Transform method of Section In all cases this requires extra memory and computational expense but it enables any equations of the general from of Equation 8 3 to be solved 2 Interaction picture Here the field components are also evolved in Fourier space using the interaction picture technique as detailed in Section 6 4 8 The form of PDE that may be solved using this technique is shown in Equation 8 5 o 3 5 x Ft k i Z a x N x a x b x x 8 4 p x H x a x b x 8 5 This form is significantly more restrictive than that of Equation 8 3 In particular the linear operator matrix must be diagonal and may not have coefficients with spatial dependence Within both of these pictures either the semi Implicit or Runge Kutta algorithms may be employed For more detail refer Section This gives rise to the matrix of algorithms listed in Table There are two main things that xmds can do to a field forward evolve integrate propagate it according to the set of PDEs and reshape filter it according to a set of functions These two actions may serve as the building blocks for an elaborate sequence of operations to be performed on the field This is done by defining a sequence of segments a segment being either an in
263. tor b was already initialised at the first lattice point in z 6 5 Discretisation and Sampling Errors In order to produce an accurate result any numerical integration must proceed with a suff ciently small time step h The difference between the correct result and the numerical result with a particular time step is known as the discretisation error For deterministic equations it is usually see next paragraph sufficient to reduce h until the result converges to the desired accuracy For stochastic equations this is not precisely correct as a greater number of time steps will induce a different choice noise from the random number generators This means that the two trajectories will be different members of the ensemble of possible trajectories so they should not be expected to converge to each other In order to verify the convergence of a stochastic integration it is necessary to ensure that the same trajectory or path is examined with time steps of differing size This can be done by averaging the noise contributions of a fine grained path to generate the noise contributions to a coarse grained path Many problems represent some localised interaction in a psuedo infite space In other words any radiation shed from the interaction should propagate away without ever being 6 5 DISCRETISATION AND SAMPLING ERRORS 113 reflected back into the region of interest In such cases boundaries must be sufficiently far away from the importan
264. tors emacs vim e Linux Distributions Gentoo Linux Redhat Linux e Other Unix and Unix like environments True64 CygWin e Scripting tools and languages aap perl python e Documentation tools doxygen TEX latex2html e Libraries fftw e Organisations Sourceforge net APAC Feedback Yes we want feedback If you have any comments about xmds and or this manual such as inaccuracies possible improvements new features what it does well etc then please email one of the current developer or the xmds web page webmaster You can find the addresses of both of these people on the xmds web page http www xmds org And please feel free to mention anything no matter how small It would be great to see xmds improve the way people want and for it to be documented the way the xmds user community wants vi Contents xiii XV 1 1 Starting from scratch 3 in m 3 1 1 1 The simple beginnings of a simulation 3 XS Bee Mi M MERE ee Be ee re 5 1 1 3 The feld l ments e qe qoe eh ORE OR 2262442 84 7 1 1 4 The sequence element 10 1 1 5 The integrate clement cont a 4 GR x RN 11 1 1 6 The output element ia x mk RR he 024 de AO sax 13 de 16 1 1 7 1 Matlab and Octave ehm Rog e Sek Ass ee bal ed 17 MAT Scilab sos LA rica REIN d de ee xd 19 T 20 de ee xus ee neo 21 Tp 22 1 2 3 Describing the HE sun Lis oe 24 6 04 24 2 0 en era 22 jor esu uA E 24
265. u ca bag basa a a csay a cgay x x a4 az h ay e enh s ki q ar h F N x F ar 8 ar e aaa hes aL a a Gay a a cay ar 1 b51 c1 ar b51 c1a bs2ag b53 b5103 c1 a bsa b5104 c1 aL x 22 ag a4 h a5 a2 hL 2 k ar ay ar h F N x F ar e as a2 hC ak ar aL a a Car Note cs 0 ar fiar fax fgay faa fsar fea x a ag as h lt h 0k ar e 26 22 P Nery ar ar h FIN x FT laz E 104 NUMERICAL MODELLING THEORY 37 ar e7 ae a2 hC aki ar 38 a a ceay 39 a a cgay 40 x x ag as h 41 a ec h ki y 49 a eU a2 h x k1 J a The f coefficients are listed in section 6 4 5 Unlike its fixed step counterpart RKAIP this algorithm requires significant memory for the fastest calculation of the Fourier space propagation as the arguments of the exponentials contain different factors of a az for each application If this turns out to be a problem XMDS offers the possibility to sacrifice speed for less memory consumption by the use of an appropriate flag As the two computed solutions are of fourth and fifth order only if the derivative operators are independent of the propagation dimension this algorithm cannot be used at all for problems where this is not the case 6 4 12 The Ninth Order Rung
266. ult_value required lt default_value gt char lt default_value gt xmds 1 2 Contains string Subelement of Path to tag lt simulation gt lt argv gt lt arg gt lt default value gt Description The default value of the command line argument Values will be converted from the string entered into the type given in the lt type gt tag If the variable is not specified at the command line then this is the value used by the simulation Example 190 LANGUAGE REFERENCE lt simulation gt lt argv gt lt arg gt lt name gt nconst lt name gt lt type gt int lt type gt lt default_value gt 2 lt default_value gt lt arg gt lt argv gt lt simulation gt 11 21 field required lt field gt xmds tags lt field gt Contains Subelement of Path to tag lt simulation gt lt globals gt lt field gt Description Container element to hold the tags describing the field to be integrated At present only one field is permitted in a simulation Example lt simulation gt lt field gt l semel weas SSS lt field gt lt simulation gt 11 21 1 name field optional lt name gt string lt name gt Contains string Subelement of Path to tag lt simulation gt lt field gt lt name gt Description The name of the field to integrate Defaults to main Example 11 21 FIELD 191 lt simulation gt lt field gt lt name gt main lt name gt lt field
267. ve on memory and also size of output file it is better to sample on some reduced lattice rather than at every point in the main field lattice It may also be desirable when sampling for a particular group of moments to transform the field to Fourier space in one some or all of the transverse dimensions At this point it is also possible to collapse one or more of the transverse dimensions by requesting either to sample a cross section or else to integrate over the dimension in question with a particular kernel function Once all field propagations are finished the sampled output can be post processed to transform the propagation dimension into Fourier space and again place the remaining transverse dimensions into any partial Fourier space which may be the same or different to the one in which the field was sampled in The original sampled moments may be used here or moments of these moments may be specified Finally if the evolution of the field involves stochastic terms it will usually be desired to perform this sequence of operations a number of times re initialising the field each time so as to determine averages and standard errors in those averages for the output moments over multiple trajectories or paths Diagrammatically the functionality of xmds is shown in Figure Usually using xmds is much simpler than this look at the worked examples 8 1 Installing and Running xmds This software package is designed to install onto Unix and Linux
268. vel syntax we have attempted to employ some degree of foresight the future is likely to see more interconnection of computer networks around the world and information is already taking a standardised form for transfer Until now this has primarily been HTML Hyper Text Mark up Language but now there is increasing interest in an extensible mark up language of which HTML is only a subset This is known as XML eXtensible Mark up Language 9 0 Also a standard for data interchange has been developed at Caltech 2 which is a subset of XML XSIL or eXtensible Scientific Interchange Language and this has been chosen as the format for field data input output 117 118 DEVELOPMENT AND PROGRAM STRUCTURE Using XML for the input script has another significant advantage the format for the simulation data takes on a tree like structure which can be extended and made as complex as is necessary This is the extensible property of XML Appropriately named element tags are used to mark up the essential data and may be nested accordingly to group the information The result is a logical human readable script with data represented in a manner preserving both synthetic and sequential relationships xmds has undergone a few development iterations for the purpose of refinement and extension There have been two particularly difficult aspects The first was in defining an appropriate equation syntax This syntax had to be in C code style since xmds mere
269. was 4 408207e 11 Time elapsed for simulation is 10 seconds where the simulation has taken approximately 10 seconds to complete 2 4 Wisdom The lt use_wisdom gt tag is the way to enable FFTW s wisdom feature This tag expects a boolean argument and by default is set to no However when set to yes you can expect an immense increase in the startup speed of your simulations Wisdom is the name fftw gives to stored information about their Fourier transform plans What fftw does before it decides to use a particular method for calculating the Fourier transform is to run some calculations beforehand to see which of the methods is the fastest this can be related to your system s architecture the size of the problem etc and then it can optionally store this information so that fftw doesn t have to go through all of the 36 EXTRA AND ADVANCED FEATURES hard work again and therefore make use of the stored wisdom about the problem at hand Enabling wisdom means that subsequent runs of the simulation will start up and run overall much faster xmds requires a place to save this accumulated wisdom so that it can be reloaded in subsequent simulation runs The way xmds does this is to save the widsom in a file called lt hostname gt wisdom where lt hostname gt is the name of the computer you are running the simulation on Note that for simulations using MPI that the wisdom filename uses the format lt hostname gt lt rank gt wi
270. would wish to have is their own personal expert computer programmer who takes their high level description of the problem and writes a low level program specifically dedicated to solving it Further it is essential that this process is error free Presented here is a solution to these difficulties xmds is a computer program that generates computer programs xmds interprets a high level description of a problem and in turn writes a low level computer program that a compiler can compile and optimise The executable file produced then solves the problem as quickly and efficiently as possible Further once xmds is debugged the process of code generation at least from the high level script onwards becomes error free Of course xmds cannot be written to be able to process absolutely any problem xmds has its own scope of capability which is outlined in Chapter We begin by outlining some of the basic theory involved in numerical modelling and then follow on with the develpoment strategy behind xmds and the structure of the source code We then look at the functionality of xmds and illustrate this with worked examples We close with a brief outlook for future development Deferential equations Numerical Modelling Theory 6 1 Differential Equations The majority of numerical models for real world problems are based on one or more differ ential equations involving the various parameters concerned The large majority of which may be written
271. x algorithms that xmds currently contains are RK4EX Fourth order Runge Kutta in the explicit picture RK4IP Fourth order Runge Kutta in the interaction picture ARK45EX adaptive step size Runge Kutta Fehlberg in the explicit picture ARK45IP adaptive step size Runge Kutta Fehlberg in the interaction picture SIEX Semi implicit method in the explicit picture SIIP Semi implicit method in the interaction picture Example lt simulation gt lt sequence gt lt integrate gt lt algorithm gt RK4IP lt algorithm gt lt integrate gt lt sequence gt lt simulation gt 11 22 2 2 interval required lt interval gt double lt interval gt Contains double Subelement of Path to tag lt simulation gt lt sequence gt lt integrate gt lt interval gt Description The integration range This is the interval over which the main propagation dimension will be integrated Example lt simulation gt lt sequence gt lt integrate gt lt interval gt 20 lt interval gt lt integrate gt lt sequence gt lt simulation gt 200 LANGUAGE REFERENCE 11 22 2 3 iterations optional lt iterations gt int lt iterations gt Contains integer Subelement of Path to tag lt simulation gt lt sequence gt lt integrate gt lt iterations gt Description When using one of the semi implicit algorithms this option can be altered to control the number of iterations of the algorithm
272. xample lt simulation gt lt sequence gt lt integrate gt lt smallmemory gt yes lt smallmemory gt lt integrate gt lt sequence gt lt simulation gt 11 22 2 8 cutoff optional lt cutoff gt double lt cutoff gt Contains double Subelement of Path to tag lt simulation gt lt sequence gt lt integrate gt lt cutoff gt Description If one of the adaptive step size methods is used the cutoff value sets the threshold for the function values that are included in the determination of relative errors Grid points where the function is less than cutoff peakvalue are not included The value defaults to 1 1000 Example 11 22 SEQUENCE 203 lt simulation gt lt sequence gt lt integrate gt lt cutoff gt 1e 5 lt cutoff gt lt integrate gt lt sequence gt lt simulation gt 11 22 2 9 halt non finite optional halt non finite bool lt halt non finite Contains boolean Subelement of Path to tag simulation sequence integrate halt non finite Description If yes then halt the current integration pass if it produces a non finite num ber in the first component of the main vector Since non finite numbers usually pro pogate quickly through the components of the main vector this will quickly detect any non finite number Examples of non finite numbers are 1 0 0 0 which is infinite and 0 0 0 0 which is not a number NaN Some platforms proc

xmds: e X tensible M ulti D imensional S imulator

Contents

Download Pdf Manuals

Related Search

Related Contents