Hi David,
thanks for the info. Anyway, I would recommend you not to use a single named integral for both the volume and surface integration, as I have written below, at least for the moment.
r.
On 08/30/2012 10:01 AM, David Libault wrote:
Hi Robert,
No I did not get this error, and I was going to tell you... My platform is Mac OSX + MacPorts. I use MacPorts python 2.7 and sfepy installed manually from git.
Le mercredi 29 août 2012 18:07:50 UTC+2, Robert Cimrman a écrit :
On 08/29/2012 05:59 PM, Robert Cimrman wrote:
Hi David,
On 08/29/2012 05:21 PM, David Libault wrote:
Hi Robert,
Enclosed are a mesh file and a problem file for the 2D electrostatic problem without volume charge of a square of 1x1 on which 1 V is applied on 2 opposite sides (Gamma_Left and Gamma_Right), and the default 0 Von Neumann boundary conditions on the 2 other sides.
The Electric field is evaluated as the gradient of V, and looks good (Ex is small, and Ey = 1 over the entire domain).
Yes, it looks good.
Flux 1 and Flux 2 evaluate the flux intensity of current with conductivity sigma=1 across Gamma_Left and Gamma_Right. They should have opposite sign and absolute value of 1.
It looks that the d_surface_flux is computed on each element of the surface so the returned value is an array, but the values in the array are all 1... Shouldn't they be 1/(element surface) so that the sum of the array gives the total flux over the all Gamma_Left or Gamma_Right surface ?
There is a glitch in using integrals I noticed when preparing the example for Alec: once you use an integral for a volume integration (volume term), you should not use it for a surface term. If you do that, you should get an error - and indeed, there is:
fmf_mulAB_nn(): ERR_BadMatch: (2 2 1) == (3 2 2) * (2 2 1)
Ha, you actually did not get this error, right? What's your platform and Python version?
printed just before your flux output. This should lead to raising an exception, but for some reason it does not, the exception is raised only when I step over the relevant code piece with a debugger. If I just run the code, it silently continues, and leads to wrong results you got. So this is definitely a bug, but I am not sure why it behaves this way. It used to work.
So to fix you problem: replace flux1 = pb.evaluate('d_surface_flux.i1.Gamma_Left(coef.sigma, v)', mode='el_avg')
with
flux1 = pb.evaluate('d_surface_flux.2.Gamma_Left(coef.sigma, v)', mode='el_avg')
and do the same with flux2. The number is the quadrature order... Then you should get the "correct" fluxes. If fact, there is another problem - 'el_avg' mode divides the resulting integrated value in an element by the volume/area/length of that element. That's why you get 1 or -1, and the sums 9 and -9, as there are 9 edges on each of Gamma_Left, Gamma_Right.
So it seems another evaluation mode is needed (let's say called 'el') - one that just integrates a value over an element, without computing the average by dividing by volume. Do you agree?