Toon put in the basic generator (fcode) and Oyvind worked on doing the
fuller generator but I don't think it was finished.  We should
probably have a F95CodeGen and F77CodeGen to separate out the major
differences, but I'm not a Fortran guy so I didn't pay attention to
the minutia when it was built.

-- Andy



2011/11/8 Ondřej Čertík <[email protected]>:
> Hi,
>
> I finally needed to generate Fortran code from a SymPy expression. I
> couldn't find any documentation for it. Is there any? I found the
> FCodeGen object, but
> it has no useful docstring.... So I looked at tests and those are
> good. Here is my real life code that I needed today:
>
>
> from StringIO import StringIO
> from sympy.utilities.codegen import FCodeGen, Routine
> from sympy import var, legendre, integrate, exp
> var("l r12 r1 r2 D")
> f = (2*l+1) / (2*r1*r2) * integrate(exp(-r12/D)*legendre(l,
> (r1**2-r12**2+r2**2) / (2*r1*r2)), (r12, r1-r2, r1+r2))
> routines = []
> for _l in range(3):
>    expr = f.subs(l, _l).doit().simplify()
>    routines.append(Routine("V%d" % _l, expr))
> code_gen = FCodeGen()
> output = StringIO()
> code_gen.dump_f95(routines, output, "file", False, False)
> source = output.getvalue()
> print source
>
>
> This prints:
>
>
> REAL*8 function V0(D, r1, r2)
> implicit none
> REAL*8, intent(in) :: D
> REAL*8, intent(in) :: r1
> REAL*8, intent(in) :: r2
> V0 = (-D*exp(-(r1 + r2)/D) + D*exp(-(r1 - r2)/D))/(2*r1*r2)
> end function
> REAL*8 function V1(D, r1, r2)
> implicit none
> REAL*8, intent(in) :: D
> REAL*8, intent(in) :: r1
> REAL*8, intent(in) :: r2
> V1 = 3*D*(-D**2*exp(2*r2/D) + D**2 - D*r1*exp(2*r2/D) + D*r1 + D*r2*exp( &
>      2*r2/D) + D*r2 + r1*r2*exp(2*r2/D) + r1*r2)*exp(-r1/D - r2/D)/(2* &
>      r1**2*r2**2)
> end function
> REAL*8 function V2(D, r1, r2)
> implicit none
> REAL*8, intent(in) :: D
> REAL*8, intent(in) :: r1
> REAL*8, intent(in) :: r2
> V2 = 5*D*(9*D**4*exp(2*r2/D) - 9*D**4 + 9*D**3*r1*exp(2*r2/D) - 9*D**3* &
>      r1 - 9*D**3*r2*exp(2*r2/D) - 9*D**3*r2 + 3*D**2*r1**2*exp(2*r2/D &
>      ) - 3*D**2*r1**2 - 9*D**2*r1*r2*exp(2*r2/D) - 9*D**2*r1*r2 + 3*D &
>      **2*r2**2*exp(2*r2/D) - 3*D**2*r2**2 - 3*D*r1**2*r2*exp(2*r2/D) - &
>      3*D*r1**2*r2 + 3*D*r1*r2**2*exp(2*r2/D) - 3*D*r1*r2**2 + r1**2*r2 &
>      **2*exp(2*r2/D) - r1**2*r2**2)*exp(-r1/D - r2/D)/(2*r1**3*r2**3)
> end function
>
>
>
>
> It looks pretty good. I would like to use "real(dp)" instead of
> "REAL*8", as well as probably create the whole module and just put
> "implicit none" there just one time. But otherwise it looks good, it
> compiles and does the job. Nice!
>
> I just found the documentation here:
>
> http://docs.sympy.org/dev/modules/utilities/codegen.html
>
> So I think I should use the codegen() function. Here is how to use it
> to do the same as above:
>
>
> from sympy.utilities.codegen import codegen
> from sympy import var, legendre, integrate, exp
> var("l r12 r1 r2 D")
> f = (2*l+1) / (2*r1*r2) * integrate(exp(-r12/D)*legendre(l,
> (r1**2-r12**2+r2**2) / (2*r1*r2)), (r12, r1-r2, r1+r2))
> routines = []
> for _l in range(3):
>    expr = f.subs(l, _l).doit().simplify()
>    routines.append(("V%d" % _l, expr))
> [(f_name, f_code), (interface_name, f_interface)] = codegen(routines,
> "F95", "potential", header=False)
> print f_code
>
>
> The f_interface is the f95 interface to the functions above. I don't
> use such things in Fortran, I just create a module and then import it
> from the module using "use moduly, only: V0". See here for more
> details:
>
> https://github.com/qsnake/qsnake/wiki/Fortran-best-practices
>
> But the f_code does what it should. So I think I'll add the
> functionality in there, as an optional parameter, to create a regular
> Fortran module.
>
> Ondrej
>
> --
> You received this message because you are subscribed to the Google Groups 
> "sympy" group.
> To post to this group, send email to [email protected].
> To unsubscribe from this group, send email to 
> [email protected].
> For more options, visit this group at 
> http://groups.google.com/group/sympy?hl=en.
>
>

-- 
You received this message because you are subscribed to the Google Groups 
"sympy" group.
To post to this group, send email to [email protected].
To unsubscribe from this group, send email to 
[email protected].
For more options, visit this group at 
http://groups.google.com/group/sympy?hl=en.

Reply via email to