question about the whereNegative

Asked by Lou Jichao

Dear all,
How could I do, if I want to know how many gauss points in the circle region with radius r0.
What I have done is'
gx = Function(mydomain).getX()
region = r0-length(gx-ic)
nt = whereNegative(region)'
But I don't know how to print how many points where nt is equal to 0
Thank you very much in advance for your help or suggestions.
Lou Jichao

Question information

English Edit question
esys-escript Edit question
No assignee Edit question
Solved by:
Adam Ellery
Last query:
Last reply:
Revision history for this message
Lou Jichao (ljc-1) said :

And what should I do if I want to know the nt value of with the number 5 gauss point.

Revision history for this message
Best Adam Ellery (aellery) said :

If the circle is centered on (x0,y0) and has radius r0 then it is described by the equation

(x-x0)^2 + (y-y0)^2 = r0^2,

where x and y are the spatial coordinates, and all interior points satisfy the inequality

(x-x0)^2 + (y-y0)^2 - r0^2 < 0.

You can find these points using something like the following:

gx = Function(mydomain).getX()
tmp = (gx-ic)*(gx-ic)
tmpB = tmp[0]+tmp[1]-r0*r0
tmpC = whereNegative(tmpB)
answer = tmpC.toListOfTuples()

Here, tmp[0] contains the expression (x-x0)^2, tmp[1] contains the expression (y-y0)^2, tmpB contains the expression (x-x0)^2 + (y-y0)^2 - r0^2 and tmpC is a Data object that is equal to unity if the corresponding point is inside the zero and is zero otherwise. If you are debugging your code, you can check the values of these Data objects using the toListOfTuples() command (e.g. using the command tmp.toListOfTuples()).

The variable `answer' contains the same values as tmpC but as a list. The number of points inside the circle is then given by sum(answer) and if you wish to know if the 5th point, for example, is inside the circle, you can use the command answer[5].

Revision history for this message
Lou Jichao (ljc-1) said :

Thanks Adam Ellery, that solved my question.