boundary elements error

[mari sadat]

hello every one

I am using iso2mesh to mesh a 3 layer cube with a cylinder in it. The cylinder is in the first and second layer.

i wrote a code that seems very well in plotting each part but i have  a trouble when i run it for monte carlo simulating.

the error is this: Ill defined mesh - some boundary elements are trapped between two elements. A boundary element must belong to a single element. Check Mesh.

I dont know what is my mistake. Do you have any idea to help me?

%% 3layer cube

[no,fc,c0]=latticegrid([0 3],[0 3],[0 2 3.7 4]);

[no,fc]=removeisolatednode(no,fc);

fc2=cell2mat(fc);

fc=[fc2(:,[1 2 3]); fc2(:,[1 3 4])];  % convert the square facets to triangles

%figure; plotmesh(no,fc);

%mesh parameter of each layer

c0(:,4)=[0.001;0.001;0.001];

%mesh the 3layer cube and obtaining node & face & element

[node,elem,face]=surf2mesh(no,fc,[],[],1,[],c0);

[node,elem]=removeisolatednode(node,elem);

%figure; plotmesh(node,elem);

%% add the cylinder to cube

c1=[1.5 1.5 3.5];

c2=[1.5 1.5 3.9]; 

[cnode,cface,celem]=meshacylinder(c1,c2,.2);

%[cnode,cface]=meshcheckrepair(cnode,cface,'dup');

%merding the cube and cylinder elements

[newnode,newelem]=mergemesh(node,elem,cnode,celem);

%figure;plotmesh(newnode,newelem,'y>1.5');

%figure;plotmesh(newnode,face,'y>1');

[Qianqian Fang]

On 10/20/21 12:42 PM, mari sadat wrote:

hello every one

I am using iso2mesh to mesh a 3 layer cube with a cylinder in it. The cylinder is in the first and second layer.

i wrote a code that seems very well in plotting each part but i have  a trouble when i run it for monte carlo simulating.

the error is this: Ill defined mesh - some boundary elements are trapped between two elements. A boundary element must belong to a single element. Check Mesh.

I dont know what is my mistake. Do you have any idea to help me?

hi Mari,

mergemesh can not resolve/handle overlapping meshes. it only concatenates spatially disconnected surface or tetrahedral meshes without checking whether they intersect or not.

if you are merging surface meshes, you should use surfboolean or mergesurf instead to resolve intersections.

https://groups.google.com/g/iso2mesh-users/search?q=mergemesh

try adding these lines at the bottom of your script and see the result (attached)

[no1,fc1]=mergesurf(no,fc,cnode,cface);
[newnode,newelem]=s2m(no1,fc1,1,10,'tetgen1.5',c0);
plotmesh(newnode,newelem,'y>1.5')

you can adjust your seeds in c0 to label regions accordingly.

Qianqian

%% 3layer cube

[no,fc,c0]=latticegrid([0 3],[0 3],[0 2 3.7 4]);

[no,fc]=removeisolatednode(no,fc);

fc2=cell2mat(fc);

fc=[fc2(:,[1 2 3]); fc2(:,[1 3 4])];  % convert the square facets to triangles

%figure; plotmesh(no,fc);

%mesh parameter of each layer

c0(:,4)=[0.001;0.001;0.001];

%mesh the 3layer cube and obtaining node & face & element

[node,elem,face]=surf2mesh(no,fc,[],[],1,[],c0);

[node,elem]=removeisolatednode(node,elem);

%figure; plotmesh(node,elem);

%% add the cylinder to cube

c1=[1.5 1.5 3.5];

c2=[1.5 1.5 3.9]; 

[cnode,cface,celem]=meshacylinder(c1,c2,.2);

%[cnode,cface]=meshcheckrepair(cnode,cface,'dup');

%merding the cube and cylinder elements

[newnode,newelem]=mergemesh(node,elem,cnode,celem);

%figure;plotmesh(newnode,newelem,'y>1.5');

%figure;plotmesh(newnode,face,'y>1');

--
You received this message because you are subscribed to the Google Groups "iso2mesh-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to iso2mesh-user...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/iso2mesh-users/bdd6e2ef-770c-4796-a579-ef6124e5672en%40googlegroups.com.

[mari sadat]

hi dear dr.Fang

First of all i appreciate you for such a nice meshing toolbox  and also thank you very much for your fast reply.

your solution about mergesurf function is great. Ididn't familiar with this function.

i did so but unfortunatly the same error about boundary appeared.

i tried different approches but non of them reply.i don't know what i should do.

thanks dr.Fang

best regards

[Qianqian Fang]

On 10/21/21 3:18 AM, mari sadat wrote:

hi dear dr.Fang

First of all i appreciate you for such a nice meshing toolbox  and also thank you very much for your fast reply.

your solution about mergesurf function is great. Ididn't familiar with this function.

i did so but unfortunatly the same error about boundary appeared.

i tried different approches but non of them reply.i don't know what i should do.

make sure you assigned the new mesh to mmc simulation. you may want to remove your mergemesh call, and do a clear all.

if the error still persists, please send me the full mmc simulation script.

To view this discussion on the web visit https://groups.google.com/d/msgid/iso2mesh-users/ba406a58-9de1-42bd-a007-b3e9ef45c585n%40googlegroups.com.

[mari sadat]

i checked it many times but i get the same error. i'm using ValoMc for the monte carlo simulating not mmc.

this is my full script.

Thanks for your attention dr.Fang.

clc

close all

clear all

%% 3layer cube

[no,fc,c0]=latticegrid([0 3],[0 3],[0 2 3.7 4]);

fc2=cell2mat(fc);

fc=[fc2(:,[1 2 3]); fc2(:,[1 3 4])];  % convert the square facets to triangles

%figure; plotmesh(no,fc);

%mesh parameter of each layer

c0(:,4)=[0.001;0.001;0.001];

%% add the cylinder to cube

c1=[1.5 1.5 3.5];

c2=[1.5 1.5 3.9]; 

[cnode,cface,celem]=meshacylinder(c1,c2,.2);

[cnode,cface]=meshcheckrepair(cnode,cface,'dup');

[no1,fc1]=mergesurf(no,fc,cnode,cface);

%figure;plotmesh(no1,fc1,'y>1')

[newnode,newelem,newface]=s2m(no1,fc1,1,10,'tetgen1.5',c0);

%figure;plotmesh(newnode,newelem,'y>1.5');

%figure;plotmesh(newnode,newface,'y>1');

vmcmesh.r=newnode;

vmcmesh.BH=newface(:,1:3);

vmcmesh.H=newelem(:,1:4);

%% accessing the cylinder & layers elements and display them

layer1=find(newelem(:,5)==1);    

layer2=find(newelem(:,5)==2);

layer3=find(newelem(:,5)==3);

cylindera=find(newelem(:,5)==4);

cylinderb=find(newelem(:,5)==5);

%% Create a medium and a boundary for the mesh

vmcmedium.scattering_anisotropy=0.9;

vmcmedium.refractive_index=1.37;

vmcmedium.absorption_coefficient(layer1)=0.03;

vmcmedium.scattering_coefficient(layer1)=17.3;

vmcmedium.absorption_coefficient(layer2)=0.008;

vmcmedium.scattering_coefficient(layer2)=35.7;

vmcmedium.absorption_coefficient(layer3)=0.3;

vmcmedium.scattering_coefficient(layer3)=37.6;

vmcmedium.absorption_coefficient(cylindera)=2.4;

vmcmedium.scattering_coefficient(cylindera)=31.3;

vmcmedium.absorption_coefficient(cylinderb)=2.4;

vmcmedium.scattering_coefficient(cylinderb)=31.3;

vmcmedium=createMedium(vmcmesh,vmcmedium);

vmcboundary=createBoundary(vmcmesh,vmcmedium);

vmcboundary.exterior_refractive_index = 1;

%% Find boundary elements

lightsource1 = findBoundaries(vmcmesh, 'direction', [1.5 1.5 4], [1.5 1.5 5], .5);

vmcboundary.lightsource(lightsource1) = {'direct'};

%% displayes the mesh

figure('Name','meshplot')

hold on

%displayes the surrface mesh

trimesh(vmcmesh.BH,vmcmesh.r(:,1),vmcmesh.r(:,2),vmcmesh.r(:,3),'facecolor', 'r','FaceAlpha',0.1) 

% Show the cylinder 

%tetramesh(vmcmesh.H(cylindera,:), vmcmesh.r);

%tetramesh(vmcmesh.H(cylinderb,:), vmcmesh.r);

% Highlight the location for the lightsource for the plot

trimesh(vmcmesh.BH(lightsource1,:),vmcmesh.r(:,1),vmcmesh.r(:,2),vmcmesh.r(:,3),'facecolor', 'b');

xlabel('x [mm]');

ylabel('y [mm]');

zlabel('z [mm]');

view(70,40);

hold off

%% run the simulation

options.photon_count=1e8;

solution = ValoMC(vmcmesh, vmcmedium, vmcboundary ,options);

%% Visualize the solution

figure('Name','solution plot')

hold on

halfspace_elements = findElements(vmcmesh, 'halfspace', [0 0 0], [0 1 0]);

tetramesh(vmcmesh.H(halfspace_elements,:), vmcmesh.r, solution.element_fluence(halfspace_elements));

view(-10,10);

xlabel('x [mm]');

ylabel('y [mm]');

zlabel('z [mm]');

c = colorbar;

c.Label.String = 'Fluence [(W/mm^2)]';

[Qianqian Fang]

On 10/21/21 6:09 PM, mari sadat wrote:

i checked it many times but i get the same error. i'm using ValoMc for the monte carlo simulating not mmc.

this is my full script.

it appears that the error you saw came from ValoMC

https://github.com/InverseLight/ValoMC/blob/5f5928f202196f36d62b1376dd931f60ffad0a7b/cpp/Errors.hpp#L51

if s2m successfully returned a mesh and you were able to plot it, I would say the issue is not on iso2mesh's side.

You should reach out to the authors of ValoMC and see if you had properly converted your input. unfortunately I am not familiar with that code.

Qianqian

To view this discussion on the web visit https://groups.google.com/d/msgid/iso2mesh-users/7754bc5b-ea53-43ff-a14b-9ffe1072facbn%40googlegroups.com.

[mari sadat]

Thanks a lot  dr.Fang for your attention.

Kind regard.

[mari sadat]

hi dr.fang

I change my script because I need only the elements of the cylinder not its surface.

[no,fc,c0]=latticegrid([0 3],[0 3],[0 2 3.7 4]);

fc2=cell2mat(fc);

fc=[fc2(:,[1 2 3]); fc2(:,[1 3 4])];  % convert the square facets to triangles

%figure; plotmesh(no,fc);

%mesh parameter of each layer

c0(:,4)=[0.001;0.001;0.001];

%mesh the 3layer cube and obtaining node & face & element

[node,elem,face]=s2m(no,fc,1,10,'tetgen1.5',c0);

%figure; plotmesh(node,elem)

%% add the cylinder to cube

c1=[1.5 1.5 3.5];

c2=[1.5 1.5 3.9];

[cnode,cface,celem]=meshacylinder(c1,c2,.2);

[cnode,cface]=meshcheckrepair(cnode,cface,'dup');

%merding the cube and cylinder elements

[newnode,newelem]=mergemesh(node,elem,cnode,celem);

%figure;plotmesh(newnode,newelem,'y>1');

    

vmcmesh.r=newnode;

vmcmesh.BH=face(:,1:3);

vmcmesh.H=newelem(:,1:4);

I chech my script with the  author of ValoMC  and understand that I should change the vmcmesh.BH.

I let face=vmcmesh.BH.  The face include the cube walls plus horizontal floors inside the cube that separate the layers. I need only the cube walls  as BH not the floors in it.

ValoMC author said me only the cube walls should be in BH, and nothing that is inside the cube.

I should access the cube walls and I don't know how I can do this.

Do you have any idea to help me?

  

thanks alot

best regards

[Qianqian Fang]

On 10/23/21 3:41 PM, mari sadat wrote:

...

%merding the cube and cylinder elements

[newnode,newelem]=mergemesh(node,elem,cnode,celem);

%figure;plotmesh(newnode,newelem,'y>1');

    

vmcmesh.r=newnode;

vmcmesh.BH=face(:,1:3);

vmcmesh.H=newelem(:,1:4);

I chech my script with the  author of ValoMC  and understand that I should change the vmcmesh.BH.

I let face=vmcmesh.BH.  The face include the cube walls plus horizontal floors inside the cube that separate the layers. I need only the cube walls  as BH not the floors in it.

ValoMC author said me only the cube walls should be in BH, and nothing that is inside the cube.

as I mentioned in my previous reply, you SHOULD NOT use mergemesh!. instead, please use mergesurf as in the example I provided. otherwise your mesh is invalid.

also, I took a very quick look at the data structures used in ValoMC, if you need the exterior boundary surface of the mesh, you should call volface(newelem(:,1:4))

Qianqian

To view this discussion on the web visit https://groups.google.com/d/msgid/iso2mesh-users/7fb997db-2c9d-4493-8cdf-67678d664babn%40googlegroups.com.

[mari sadat]

Dear professor, I apologize to you.

Yes, you told me, but I thought because I didn't need the cylinder surfaces,  I could use mergemesh!

I correct my script with mergsurf  and use volface to access BH.

Every thing is ok and my simulation is started.

Thank you very much for your nice toolbox and also for your great help.

sincerely yours

kind regard

mari

[mari sadat]

Hi professor fang.

Excuse me for my too many questions🙏

My simulation is done and valomc is finished in about 8 hours. But plot of solution is still runnig(about11hours) and matlab is busy and i can not see solution plot yet.

Is there a problem with it?

I use  tetramesh for plot the solution.

My total elements are about 452000.

Before running montecarlo when i wont to display layers with tetramesh it tooks long time and i paused it but in order to understand the fifth column of elements and distinguish layers for adjusting optical parameters i choose 1000 elements of each of them and plot with tetramesh.

You received this message because you are subscribed to a topic in the Google Groups "iso2mesh-users" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/iso2mesh-users/Q0FAN8hsxXo/unsubscribe.
To unsubscribe from this group and all its topics, send an email to iso2mesh-users+unsubscribe@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/iso2mesh-users/e69b1491-85df-4203-ab1b-8c589b85bea4n%40googlegroups.com.

[Qianqian Fang]

On 10/26/21 1:25 PM, mari sadat wrote:

Hi professor fang.

Excuse me for my too many questions🙏

My simulation is done and valomc is finished in about 8 hours. But plot of solution is still runnig(about11hours) and matlab is busy and i can not see solution plot yet.

Is there a problem with it?

I don't know what function you used for the plotting, but if you used iso2mesh's plotmesh, this sounds too long to me.

again, if you used valomc's tools, please send your inquiries to the authors instead to get their feedback.

I use  tetramesh for plot the solution.

My total elements are about 452000.

Before running montecarlo when i wont to display layers with tetramesh it tooks long time

if you use the 3-line meshing/plotting function I sent previously (see the quoted reply I sent on Oct 20), the plotting completed in a fraction of a second.

Qianqian

To view this discussion on the web visit https://groups.google.com/d/msgid/iso2mesh-users/CAM5CT9LM8N930PDMJW7A%3DCzOxVxeL0%3DdSwSaD-iUDZgC06TViQ%40mail.gmail.com.

[mari sadat]

Hi Dr.Fang

Thanks a lot for your usefull advices.

I can plot the mesh by the plotmesh function quickly but my problem is plotting the solution.

I used tetramesh function for plotting this but it didn't respond. ValoMC authors said that Visualisation of large tetrahedral meshes is inefficient using tetramesh.

Also by using plotmesh function to visualise the solution, the matlab didn't respond. I used plotmesh in my script like this:

vmcmesh.r=newnode;

vmcmesh.BH=openedge;

vmcmesh.H=newelem(:,1:4);

plotmesh(vmcmesh.r,vmcmesh.H, solution.element_fluence,'y>1.5')

of course whenever i choose some elements boath tetramesh and plotmesh respond:

tetramesh(vmcmesh.H(1:10000,:), vmcmesh.r, solution.element_fluence(1:10000));

or

figure;plotmesh(vmcmesh.r,vmcmesh.H(1:100,1:4), solution.element_fluence(1:100,:))

for both of above lines plot is appeared but for the second this error is appeared also:

Index in position 2 exceeds array bounds (must not exceed 1).

Error in plotmesh (line 192)

   cent=meshcentroid(node,elem(:,1:4));

I don't  know how I can plot the whole solution.

Bests

Mari

 

To unsubscribe from this group and all its topics, send an email to iso2mesh-user...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/iso2mesh-users/e69b1491-85df-4203-ab1b-8c589b85bea4n%40googlegroups.com.

--
You received this message because you are subscribed to the Google Groups "iso2mesh-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to iso2mesh-user...@googlegroups.com.

[Qianqian Fang]

On 10/28/21 6:21 PM, mari sadat wrote:

Hi Dr.Fang

Thanks a lot for your usefull advices.

I can plot the mesh by the plotmesh function quickly but my problem is plotting the solution.

I used tetramesh function for plotting this but it didn't respond. ValoMC authors said that Visualisation of large tetrahedral meshes is inefficient using tetramesh.

Also by using plotmesh function to visualise the solution, the matlab didn't respond. I used plotmesh in my script like this:

vmcmesh.r=newnode;

vmcmesh.BH=openedge;

vmcmesh.H=newelem(:,1:4);

plotmesh(vmcmesh.r,vmcmesh.H, solution.element_fluence,'y>1.5')

if you were not using the correct syntax of this function, you should not expect it to work.

please read the help info carefully and use the correct parameters.

also, there are many examples in the sample/ folder if you search for this function name.

To view this discussion on the web visit https://groups.google.com/d/msgid/iso2mesh-users/e1f99330-bfe1-4757-a5ae-2a07c3b37cf9n%40googlegroups.com.

[mari sadat]

As I said before I don't have any problem with plotting the mesh. I visualise it quickly. my problem is in display the solution. 

Even if I use plotmesh(newnode,vmcmesh.H)  the matlab doesn't resppond because i omit fifth column of newelem for vmcmesh.H.

On the other hand, the dimension of my solution matrix is (number of elements*1).I think plotmesh is not suitable for my purpose.

Thanks alot for your attention

[Qianqian Fang]

On 10/29/21 2:34 PM, mari sadat wrote:

As I said before I don't have any problem with plotting the mesh. I visualise it quickly. my problem is in display the solution. 

Even if I use plotmesh(newnode,vmcmesh.H)  the matlab doesn't resppond because i omit fifth column of newelem for vmcmesh.H.

On the other hand, the dimension of my solution matrix is (number of elements*1).I think plotmesh is not suitable for my purpose.

Thanks alot for your attention

the correct plotmesh syntax for plotting element-associated values can be found in the help info

https://github.com/fangq/iso2mesh/blob/master/plotmesh.m#L19-L21

basically, you should have used

potmesh(newnode, [vmcmesh.H(:,1:4), log10(solution.elem_fluence)], 'y>1.5')

to display the results.

similarly, if the results are node-based, you should use

potmesh([newnode, log10(solution.node_fluence)], vmcmesh.H, 'y>1.5')

instead, as described here

https://github.com/fangq/iso2mesh/blob/master/plotmesh.m#L10-L11

for some reason, if your call produces extremely slow rendering and garbled elements, you should explicitly set the face-input as [], so plotmesh won't mistaken your element input as face, something like

potmesh(newnode, [], [vmcmesh.H(:,1:4), log10(solution.elem_fluence)], 'y>1.5')

To view this discussion on the web visit https://groups.google.com/d/msgid/iso2mesh-users/e32cd5cb-ab60-4b65-9b15-28e1f4a5f444n%40googlegroups.com.

[mari sadat]

Hi profssor Fang

Thank you very much for your valuable advices. I was finally able to see the answer to my simulation in a second.

I will never forget your kindness and wish you good health with happiness.

Best regards

Mari