Re: Re: Numerical evaluation is Mathematica bottleneck?!
- To: mathgroup at smc.vnet.net
- Subject: [mg69107] Re: Re: Numerical evaluation is Mathematica bottleneck?!
- From: "Jens-Peer Kuska" <kuska at informatik.uni-leipzig.de>
- Date: Wed, 30 Aug 2006 06:33:17 -0400 (EDT)
- Organization: Uni Leipzig
- References: <ed0sb8$src$1@smc.vnet.net>
- Sender: owner-wri-mathgroup at wolfram.com
Hi, and there is only the unreadable form below ?? usual I prefer equations in the traditional mathematical way -- because from you original question it is clear that one should use NDSolve[] directly, or build a SparseMatrix[] object and not the discretisation you give below. Regards Jens "gregorc" <gregor.cernivec at fe.uni-lj.si> schrieb im Newsbeitrag news:ed0sb8$src$1 at smc.vnet.net... | Hello, | | The problem consists of three coupled nonlinear equations: the Poisson | equation and the countiniuty equations for electrons and holes. The | countiniuty equations consist of convection-diffusion type equations for | electrons and hole current and the generation-recombination terms. | | The Poisson equation for non-uniform mesh at node 2: | | FV[2]=(2*(-((-V[[1]] + V[[2]])/h[[1]]) + (-V[[2]] + | V[[3]])/h[[2]]))/(h[[1]] + h[[2]])-(q*(nt[2] + n[[2]] + NA[[2]] - | ND[[2]] - p[[2]]-pt[2]))/eps[[2]] | | The electron and hole currents between nodes 1-2 and 2-3 and the | countinuity equations, respectively: | | Jn[1]=(Dn[[1]]*(-(B[(V[1] - V[2])/UT[[1]]]*n[[1]]) + B[(-V[1] + | V[2])/UT[[1]]]*n[[2]]))/h[[1]] | Jn[2]=(Dn[[2]]*(-(B[(V[2] - V[3])/UT[[2]]]*n[[2]]) + B[(-V[2] + | V[3])/UT[[2]]]*n[[3]]))/h[[2]] | FJn[2]=-2/(h[[2]]+h[[1]])*(Jn[1]-Jn[2])-R+GL | | Jh[1]=(Dp[[1]]*(B[(-V[1] + V[2])/UT[[1]]]*p[[1]] - B[(V[1] - | V[2])/UT[[1]]]*p[[2]]))/h[[1]] | Jh[2]=(Dp[[2]]*(B[(-V[2] + V[3])/UT[[2]]]*p[[2]] - B[(V[2] - | V[3])/UT[[2]]]*p[[3]]))/h[[2]] | FJh[2]=2/(h[[2]]+h[[1]])*(Jp[1]-Jp[2])-R+GL | | ,where V[[2]] is the calculation variable for the electrical potential | at node 2, n[[2]] and p[[2]] are the calculation varibles for electron | and hole concentrations, respectively. The rest is material data, built | in and/or trapped charge, recombination-generation terms, which might be | also set to zero (R,GL) and constants. B is the Bernoulli function, | which has to be approximated with Taylor series to avoid singularity | with argument 0 - zero field (V[[1]]=V[[2]]). h[[1]], h[[2]]... is the | discretization step. | | Functions:FJn[2], FV[2], FJh[2] | variables:X={n[[1]],V[[1]],p[[1]],n[[2]],V[[2]],p[[2]],n[[3]],V[[3]],p[[3]]} | | The Jacobian: | | JFV=Map[D[FV[2],#]&,X] | JFJn=Map[D[FJn[2],#]&,X] | JFJp=Map[D[FJp[2],#]&,X] | | Functions and the coresponding Jacobians are compiled: | | JFVnum=Compile[{{n,_Real,1},{V,_Real,1},{p,_Real,1},{consts},...},Evaluate[JFV]]. | | Numerical evaluation: | | JFVvar=MapThread[JFVnum,{var`n,var`V,var`p,...consts, mat. data}] | | ,where var`n=Partition[n,3,1] and same for the V and p and the rest of | the material data and constants. | | The numerical function is Threaded over the overlapping triplets of | calculation variables and data. When the numerical function and the | coresponding jacobians are calculated, the sparse array is formed with | the SparseArray function and the system is solved with the LinearSolve. | The majority of the computational time takes the MapThread of the | compiled function on the partitioned data and varibles. The functions | and jacobians are compiled for sure, since I have checked it with | compiled_function[[-3]], which returns all real numbers. | | Regards, | | gregor. | | Jens-Peer Kuska wrote: | | >Hi, | > | >the C/C++ version may be faster -- or not, it depend on | >your programming skills and the numerical methods you are | >using. | >But as long as you don't tell more about the problem and | >how you have solved it in Mathematica, there is not defined answer. | > | >Regards | > Jens | > | >gregorc wrote: | > | > | >>Dear mathgroup members, | >> | >>I am solving a non-linear semiconductor PDE problem with finite | >>difference discretization. The Compiled algebraic equation and its | >>Jacobian take as input central point and two adjacent neighbors and | >>they are Threaded over the whole mesh. A LinearSolve and some damping | >>schemes are used to obtain next step towards the solution. Typical mesh | >>is about 300 points and the solution is usually found in 7 iterations, | >>which takes about 7 seconds of computational time.This is already too | >>slow, but since I plan to solve the problem in 2D this is much too | >>slow!!! The bottleneck is the numerical evaluation of the compiled | >>equation and its Jacobian. I am thinking of using MathLink and write the | >>numerical evaluation routine in C++, which would take a list of data and | >>return evaluated function and the Jacobian. | >> | >>Is this the right way? Could the link itself be a timing bottleneck? Any | >>other suggestions? | >> | >>Regards, | >> | >>Gregor | >> | >> | >> | > | > | >__________ NOD32 1.1727 (20060826) Information __________ | > | >This message was checked by NOD32 antivirus system. | >http://www.eset.com | > | > | > | > | > |