Symbolic calculations can be done in sympy and then directly exported to C/C++ code. To do this, use the ccode function in sympy. There is a tutorial for sympy in the python subsection of the Literature page.
sympy export to C/C++
# Global imports and variables
from sympy import *
x,y,z = symbols('x y z')
t,T = symbols('t T')
Umax, R, x0, y0 = symbols("Umax, R, x0, y0")
# 3D shear (single-vortex) velocity field
shearUx = sin(2*pi*y)*sin(pi*x)**2*cos(pi * t / T)
shearUy = -sin(2*pi*x)*sin(pi*y)**2*cos(pi * t / T)
r = ((x - x0)**2 + (y-y0)**2)**0.5
shearUz =Umax * (1 - r/ R)**2*cos(pi * t / T)
# 3D shear velocity gradient.
gradShear3DU = [[diff(shearUx,x), diff(shearUx,y), diff(shearUx,z)],
[diff(shearUy,x), diff(shearUy,y), diff(shearUy,z)],
[diff(shearUz,x), diff(shearUz,y), diff(shearUz,z)]]
print(ccode(gradShear2DU))
program output:
[[2M_PIsin(M_PIx)sin(2M_PIy)cos(M_PIx)cos(M_PIt/T), -2M_PIpow(sin(M_PIy), 2)cos(2M_PIx)cos(M_PIt/T), 0],
[2M_PIpow(sin(M_PIx), 2)cos(2M_PIy)cos(M_PIt/T), -2M_PIsin(2M_PIx)sin(M_PIy)cos(M_PIy)cos(M_PIt/T), 0],
[0, 0, 0]]