Robert, thank you for your response. I've now written a function that calculates the gradient as a field variable from a scalar field variable:
def gradient(pb,region,scalar,optorder=1): from sfepy_scb.projections import project_by_component import numpy as nm
gfield = Field.from_args('gfield', nm.float64, 2, region,
approx_order=optorder - 1)
gradient = FieldVariable('gradient', 'parameter', dfield,
primary_var_name='(set-to-None)')
cfield = Field.from_args('component', nm.float64, 1, region,
approx_order=optorder - 1)
component = FieldVariable('component', 'parameter', cfield,
primary_var_name='(set-to-None)')
ev = pb.evaluate
order = 2 * (optorder - 1)
theterm= 'ev_grad.%d.%s(%s)' % (order,region.name,scalar.name)
gradient_qp=ev(theterm, mode='qp')
gradient_qp=nm.swapaxes(gradient_qp,2,3) #make the row vectors into
columns
project_by_component(gradient, gradient_qp, component, order)
return gradient
I can then use the output to evaluate line probes:
gradc= gradient(pb,omega,c) p=LineProbe([0,0],[5,0],0) gline=p(gradc)
Note that in the second line of gradient() I've imported my own version of project_by_component. That's because the original version expects a 3-component gradient and mine has only 2 components. I had to add a check for the vector dimension:
dim = nm.shape(tensor_qp)[2]
for ic in range(dim):
My next task is to submit that as an issue.
I would also like to make a similar function to calculate flux, using something like:
theterm= 'ev_diffusion_velocity.%d.%s(s.D, %s)' % (order,region.name,
scalar.name) flux_qp=ev(theterm, copy_materials=False, mode='qp')
The problem is that I don't know how to pass the material.diffusivity in as an argument. The best is could do is pass a string 's.D' because s.D does not evaluate to anything. Any suggestions?