The following runs for me:
def post_process(out, pb, state, extend = False): import numpy as nm from sfepy.fem import Integral flux = pb.evaluate('d_surface_flux.i1.Gamma(m.K, u)', mode = 'qp') flux = extend_cell_data(flux, pb.domain, 'Gamma', is_surface =
flux_field = Field('flux', np.float64, (1,),
pb.domain.regions['Omega']) flux_var = FieldVariable('flux', 'parameter', flux_field, 1, primary_var_name='(set-to-None)')
coors = nm.array([[0.25, 0.25, 0.25]], dtype=nm.float64) weights = nm.array([1.0], dtype=nm.float64) integral = Integral('ii', coors=coors, weights=weights) flux_var.data_from_qp(flux, integral) val = flux_var() print val.min(), val.max() out.update(flux_var.create_output()) return out
As i1 has 1 udrature point, I have manually create a volume integral with a single quadrature point, so that the array dimensions match. But it is not the way the correct values should be obtained :)
The correct way would be to average from surface faces to surface vertices, and set the other vertex values to zeros.
On 09/11/2012 12:16 PM, Robert Cimrman wrote:
I meant volume integral as an Integral instance=quadrature points, not the volume term :). I will try your file when I get to my laptop.
It seems that the flux shape should be (5235, 4, 1, 1) (4 quadrature
points). So you should use a volume integral in data_from_qp(), that you
use in solving the problem, as extend_cell_data() returns the flux in
the volume elements (averaged from the surface). Hmm... But the volume integral that is used in problem definition has the "dw_" prefix and cannot be used in evaluation mode. I think only "dq_" and "ev_" terms can be used in "qp" mode. I try to used "ev_grad" term but this term returns gradient (vector). It is also very nice, but I need scalar flux. :) To be more clear I attached my standalone problem definition python script and the mesh This is probably not
optimal for surface terms (a better way would be to directly average the
face data to face vertices), but it is how things are now.
Yes, yes. This is not the problem to directly average the face elements data to the face vertices, it is about several lines of code. But I want to understand more about sfepy internals. :)