mcellteam / mcell

MCell Monte Carlo Simulator of Cellular Microphysiology
Other
33 stars 10 forks source link

Overlapping/coplanar faces can cause leaks #21

Closed jczech closed 8 years ago

jczech commented 8 years ago

This issue is being created because of a problem that Bob reported while creating a tutorial for Fick's Law. I've simplified his original MDL considerably and am including it here. The root of the problem is the fact that some faces in his model are overlapping/coplanar, which results in molecules leaking out. If two faces are identical, MCell will report an error as long as the with_checks flag isn't disabled. We've known that coplanar faces can cause problems in the past and even warn against it in the QRG (http://mcell.org/documentation/qrg/mcell3_qrg.html#avoid-coincident-meshes), but we've never taken any steps to resolve it afaik, aside from the special case where they are completely overlapping (i.e. identical).

ITERATIONS = 500
TIME_STEP = 1e-06

DEFINE_MOLECULES {
  vm { DIFFUSION_CONSTANT_3D = 1e-05 }
}

box BOX {
  CORNERS = [0.01,-0.3,-0.3],[0.09,0.3,0.3]
}

plane POLYGON_LIST {
  VERTEX_LIST
  {
    // The directions (e.g. bottom, top) used below were chosen somewhat
    // arbitrarily by looking along the X-axis (from negative to positive) with
    // the Z-axis pointing up in Blender. 
    //
    // This probably won't look right if you don't use a monospaced font
    /*
       2     3
        -----
        |\  |
        | \ |
        |__\|
       1     0
    */
    [0.01,-0.3,-0.3]  // bottom right
    [0.01,0.3,-0.3]   // bottom left
    [0.01,0.3,0.3]    // top left
    [0.01,-0.3,0.3]   // top right
  }
  ELEMENT_CONNECTIONS
  {
    // The faces and vertices directly overlap. Mcell will catch this as long
    // as the "with_checks" flag isn't set to "no".
    /*[1, 2, 3]*/
    /*[3, 0, 1]*/
    // The vertices directly overlap, but not the faces. Mcell will not catch
    // this.
    [0, 1, 2]
    [2, 3, 0]
  }
}

INSTANTIATE world OBJECT
{
  box OBJECT box {}
  box_rel RELEASE_SITE
  {
    SHAPE = world.box
    MOLECULE = vm
    NUMBER_TO_RELEASE = 250
  }
  plane OBJECT plane {}
}

VIZ_OUTPUT
{
  MODE = CELLBLENDER
  FILENAME = "./viz_data/fick"
  MOLECULES
  {
    NAME_LIST {ALL_MOLECULES}
    ITERATION_NUMBERS {ALL_DATA @ ALL_ITERATIONS}
  }
}
jczech commented 8 years ago

edit: fixed a poorly worded sentence

Here's an additional example where MCell produces an error when the faces are both coplanar and not identical (this is correct behaviour):

ITERATIONS = 1000
TIME_STEP = 1e-6

DEFINE_MOLECULES
{
  sm { DIFFUSION_CONSTANT_2D = 1e-06 }
}

Triangle POLYGON_LIST
{
  VERTEX_LIST
  {
    [ -1, -1, 0 ]
    [ 1, -1, 0 ]
    [ -1, 1, 0 ]
  }
  ELEMENT_CONNECTIONS
  {
    [ 0, 1, 2 ]
  }
}

INSTANTIATE Scene OBJECT
{
  Triangle OBJECT Triangle{}
  Triangle2 OBJECT Triangle{ SCALE = [1.5, 1.5, 1.5] TRANSLATE = [0.2, 0.2, 0.0]}
  // This version of Triangle2 will work since it's slightly offset along the
  // z-axis
  /*Triangle2 OBJECT Triangle{ SCALE = [1.5, 1.5, 1.5] TRANSLATE = [0.2, 0.2, 0.05]}*/
  rel RELEASE_SITE
  {
    SHAPE = Scene.Triangle
    MOLECULE = sm'
    NUMBER_TO_RELEASE = 1000
  }
}

This is correct behaviour and is merely included to contrast with the previous case where the triangles were both overlapping and not identical, but MCell did not catch the error.

tmbartol commented 8 years ago

Right. That's the correct behavior for this case. Are you saying this is incorrect behavior?

Tom

On Wed, Oct 14, 2015 at 07:56:56AM -0700, Jacob Czech wrote:

Here's an example where MCell catches the coplanar faces (and issues an error) even though the vertices aren't coincident at all:

ITERATIONS = 1000
TIME_STEP = 1e-6

DEFINE_MOLECULES
{
  sm { DIFFUSION_CONSTANT_2D = 1e-06 }
}

Triangle POLYGON_LIST
{
  VERTEX_LIST
  {
    [ -1, -1, 0 ]
    [ 1, -1, 0 ]
    [ -1, 1, 0 ]
  }
  ELEMENT_CONNECTIONS
  {
    [ 0, 1, 2 ]
  }
}

INSTANTIATE Scene OBJECT
{
  Triangle OBJECT Triangle{}
  Triangle2 OBJECT Triangle{ SCALE = [1.5, 1.5, 1.5] TRANSLATE = [0.2, 0.2, 0.0]}
  // This version of Triangle2 will work since it's slightly offset along the
  // z-axis
  /*Triangle2 OBJECT Triangle{ SCALE = [1.5, 1.5, 1.5] TRANSLATE = [0.2, 0.2, 0.05]}*/
  rel RELEASE_SITE
  {
    SHAPE = Scene.Triangle
    MOLECULE = sm'
    NUMBER_TO_RELEASE = 1000
  }
}

Reply to this email directly or view it on GitHub: https://github.com/mcellteam/mcell/issues/21#issuecomment-148075980

jczech commented 8 years ago

No, I agree that this second example is how MCell should behave. I suppose my wording was kind of confusing. For the sake of completeness, I just wanted to include an example showing this case:

Because my original example only showed this:

So, in short, I just wanted to show that MCell sometimes does the right thing even if the overlapping triangles aren't identical.

As a side note, some sections of the QRG aren't clear whether the burden of checking is on the user or MCell, like this one: http://mcell.org/documentation/qrg/mcell3_qrg.html#avoid-coincident-meshes. Other places are more explicit about it like the beginning of this section: http://mcell.org/documentation/qrg/mcell3_qrg.html#geometrical-objects.

tmbartol commented 8 years ago

OK. I see now. Yeah, we need to clarify this consistently in the QRG.

FYI, in MCell version 2 I had ray-tracing code that did a good job of correctly handling the special case of identical overlapping triangles. However, it is generally very difficult for the user to create and guarantee that they have created perfectly identical triangles. So in MCell version 3 we decided that it is safer to disallow them all together.

Tom

On Wed, Oct 14, 2015 at 09:35:28AM -0700, Jacob Czech wrote:

No, I agree that this second example is how MCell should behave. I suppose my wording was kind of confusing. For the sake of completeness, I just wanted to include an example showing this case:

  • non-identical but overlapping: error caught (correct behaviour)

Because my original example only showed this:

  • identical and overlapping: error caught (correct behaviour)
  • non-identical but overlapping: error not caught (incorrect behaviour)

So, in short, I just wanted to show that MCell sometimes does the right thing even if the overlapping triangles aren't identical.

As a side note, some sections of the QRG aren't clear whether the burden of checking is on the user or MCell, like this one: http://mcell.org/documentation/qrg/mcell3_qrg.html#avoid-coincident-meshes. Other places are more explicit about it like the beginning of this section: http://mcell.org/documentation/qrg/mcell3_qrg.html#geometrical-objects.


Reply to this email directly or view it on GitHub: https://github.com/mcellteam/mcell/issues/21#issuecomment-148108979

haskelladdict commented 8 years ago

Even if we properly detect collisions in this case there still is the question though how to proceed. E.g., let's say molecule M hits triangle A which is coplanar with triangle B and C and the hit point is also in B and C. Let's assume that each of them has different region properties, i.e., one is reflective to M, one absorptive and one transparent. What should happen to particle M?

In general, this may even be an issue if the triangles aren't coplanar but intersect in a line which M hits exactly. This would be the case if a sampling box intersects a mesh object (without ever being coplanar). I don't remember if we do anything smart in the code to deal with this case. We probably do since I haven't noticed any leakage in this case but perhaps we were just lucky.

tmbartol commented 8 years ago

On Wed, Oct 14, 2015 at 10:52:17AM -0700, haskelladdict wrote:

Even if we properly detect collisions in this case there still is the question though how to proceed. E.g., let's say molecule M hits triangle A which is coplanar with triangle B and C and the hit point is also in B and C. Let's assume that each of them has different region properties, i.e., one is reflective to M, one absorptive and one transparent. What should happen to particle M?

In MCell 2 we gave priority to transparency first, absorption second, reflective, third, or some such. But it is weird and hard to control. I agree with disallowing co-planar faces in MCell 3 for all these reasons.

In general, this may even be an issue if the triangles aren't coplanar but intersect in a line which M hits exactly. This would be the case if a sampling box intersects a mesh object (without ever being coplanar). I don't remember if we do anything smart in the code to deal with this case. We probably do since I haven't noticed any leakage in this case but perhaps we just were lucky.

Yes, I'm pretty sure Rex put in special code to handle this case. I have not noticed any leaks in the large neuropil model which makes extensive use of transparent counting volumes that intersect with the dendritic spines.

Tom


Reply to this email directly or view it on GitHub: https://github.com/mcellteam/mcell/issues/21#issuecomment-148134709

haskelladdict commented 8 years ago

Yes, I'm pretty sure Rex put in special code to handle this case. I have not noticed any leaks in the large neuropil model which makes extensive use of transparent counting volumes that intersect with the dendritic spines.

Good to know. I'll have to take a look to see what we do in this case. Once we understand all the cases it might be good to write these up in a policy document.

tmbartol commented 8 years ago

Yes, that would be extremely helpful. :-)

Tom

On Wed, Oct 14, 2015 at 11:13:05AM -0700, haskelladdict wrote:

Yes, I'm pretty sure Rex put in special code to handle this case. I have not noticed any leaks in the large neuropil model which makes extensive use of transparent counting volumes that intersect with the dendritic spines.

Good to know. I'll have to take a look to see what we do in this case. Once we understand all the cases it might be good to write these up in a policy document.


Reply to this email directly or view it on GitHub: https://github.com/mcellteam/mcell/issues/21#issuecomment-148141758

jczech commented 8 years ago

This has been fixed in commit 4262888.