Page - 168 - in Programming for Computations – Python - A Gentle Introduction to Numerical Simulations with Python

168 5 SolvingPartialDifferentialEquations tol = 1E-12 for i in range(0, u.shape[0]): diff = abs(u_exact(x, t[i]) - u[i,:]).max() assert diff < tol, ’diff=%.16g’ % diff print ’diff=%g at t=%g’ % (diff, t[i]) WithN D4wereproducethe linearsolutionexactly. Thisbringsconfidence to the implementation,which is justwhatwe need for attacking a real physical problem next. Problemswithreusingtherhs function Therhs functionmust takeuandtasarguments,because that is requiredbythe ode_FE function. What about the variablesbeta,dx, L,x, dsdt,g, anddudx that therhs functionneeds? Theseareglobal in the solutionwehavepresented sofar.Unfortunately, thishasanundesiredsideeffect:wecannot import therhs functioninanewfile,definedudxanddsdt in thisnewfileandget the imported rhs to use these functions. The imported rhswill use the global variables, including functions, in itsownmodule. Howcanwefindsolutions tothisproblem?Technically,wemustpacktheex- tradatabeta,dx,L,x,dsdt,g, anddudxwith therhs function,which requires moreadvancedprogrammingconsideredbeyondthescopeof this text. Aclass is thesimplest constructionforpackinga function togetherwithdata, see the beginning of Chapter 7 in [13] for a detailed example on how classes can be used in such a context. Another solution in Python, and especially in computer languages supporting functional programming, is so called closures. They are also covered in Chapter 7 in thementioned reference and behave in amagicway. The thirdsolution is toallowanarbitrarysetofarguments forrhs in a list tobe transferred toode_FEand thenback torhs.AppendixH.4 in [13] explains the technicaldetails. 5.1.4 Application:HeatConductioninaRod Let us return to the case with heat conduction in a rod (5.1)–(5.4). Assume that the rod is 50cm long andmade of aluminumalloy 6082. Theˇ parameter equals =.%c/, where is the heat conduction coefficient, % is the density, and c is the heat capacity.Wecanfindpropervalues for thesephysicalquantities in thecaseof aluminumalloy 6082: % D 2:7 103kg/m3, D 200 WmK, c D 900 JKkg. This results inˇD =.%c/D8:2 10 5m2=s. Preliminarysimulationsshowthatweare close toaconstant steadystate temperatureafter1h, i.e.,T D3600s. Therhs function fromtheprevious sectioncanbe reused,only the functionss, dsdt,g, anddudxmustbechanged(seefilerod_FE.py): def dudx(t): return 0 def s(t): return 323