[MOAB-dev] Problems with remove_adjacencies

[Iulian Grindeanu]

Hello, 

----- Original Message -----


Hi Iulian, 

I'll try and give an explanation of what I am doing a bit better. Basically, I import a mesh as a VTK. When a face meets a condition for crack insertion it is marked in a tag. Then I go about re-arranging the mesh. At this point I have put all the elements into a meshset. 

When splitting the faces I duplicate the nodes I require to (not all nodes on the face are necessarily duplicated) and I change the connectivity of the tets, faces and edges. I do this in the code shown below. 

MBRange entities=side_elems; 
result = moab->get_adjacencies(side_elems,2,true,side_elems_faces,MBInterface::UNION); 
entities.merge(side_elems_faces); 
MBRange side_elems_edges; 
result = moab->get_adjacencies(side_elems,1,true,side_elems_edges,MBInterface::UNION); 
entities.merge(side_elems_edges); 
for(MBRange::iterator eit=entities.begin(); eit!=entities.end(); eit++){ 
const MBEntityHandle* conn; 
int nb_nodes; 
result = moab->get_connectivity(*eit,conn,nb_nodes);RR 
MBEntityHandle new_conn[nb_nodes]; 
for(int ii=0;ii<nb_nodes;ii++) 
{ 
if(conn[ii]==old_node) new_conn[ii] = new_node; 
else new_conn[ii] = conn[ii]; 
} 
result = moab->set_connectivity(*eit,new_conn,nb_nodes);RR //Set new connectivity 
result = moab->remove_adjacencies(*eit,&old_node,1); 
result = moab->remove_adjacencies(old_node,&*eit,1); 
} 



So, in this piece of code, you are trying to replace one old_node with the new_node. When you set_connectivity(*eit, new_conn, ...), the up-vertex adjacencies are adjusted too, you do not have to "remove adjacency". These "removes" have no effect. I assume new_node is created already, and it is in a position close to old_node (moab does not do any quality checks at this point) 
You could get the same effect of this code by: 
1) get all entities adjacent to old_node (with create flag true, if you need to) 
2) do the for loop over those entities only 
Your entities range contain all faces, edges, elements in "side_elems" range (overkill) 

<blockquote>
With this change in connectivity some adjacencies are no longer physically able to be in the mesh so this is where I try and change things by removing adjacencies. To do this I use the following code. 


</blockquote>
I don't think there is anything to remove at his point. The only thing you can remove (actually delete) here is old_node, it is free at this point. 

<blockquote>
MBRange other_side_elems = subtract(elems,side_elems); 
for(MBRange::iterator eit = other_side_elems.begin();eit!=other_side_elems.end();eit++){ 
MBRange other_side_elems_nodes; 
result = moab->get_connectivity(&*eit,1,other_side_elems_nodes); 
for(MBRange::iterator fit = side_elems_faces.begin();fit!=side_elems_faces.end();fit++){ 
MBRange face_nodes; 
result = moab->get_connectivity(&*fit,1,face_nodes); 
MBRange intersect_nodes = intersect(face_nodes,other_side_elems_nodes); 
assert(intersect_nodes.size()<=6); 
if(intersect_nodes.size()!=6){ 
result = moab->remove_adjacencies(*eit,&*fit,1);RR 
result = moab->remove_adjacencies(*fit,&*eit,1);RR 
} 
MBRange test1,test2; 
moab->get_adjacencies(&*eit,1,2,false,test1); 
moab->get_adjacencies(&*fit,1,3,false,test2); 
} 
for(MBRange::iterator fit = side_elems_edges.begin();fit!=side_elems_edges.end();fit++){ 
MBRange edge_nodes; 
result = moab->get_connectivity(&*fit,1,edge_nodes); 
MBRange intersect_nodes = intersect(edge_nodes,other_side_elems_nodes); 
assert(intersect_nodes.size()<=3); 
if(intersect_nodes.size()!=3){ 
result = moab->remove_adjacencies(*eit,&*fit,1);RR 
result = moab->remove_adjacencies(*fit,&*eit,1);RR 
} 
MBRange test1,test2; 
moab->get_adjacencies(&*eit,1,2,false,test1); 
moab->get_adjacencies(&*fit,1,3,false,test2); 
} 
} 

Obviously there are dependencies and things that I am not showing but without showing my whole code which is quite big I couldn't really do this. Hopefully this has given you some more clarity on what I am doing. 

Regards 

Graeme 

</blockquote>
So I really don't get this second part of the code. You did not duplicate yet any face or edge. 

Let's assume that you have this simple model, 2 quads, 1452, 2563 with 6 nodes 
1 2 3 
4 5 6 
You want to add a crack at 25, so you want to end up something like this 

1 2 2' 3 
4 5 6 

You first create a new node, 2'. This node is, so far, free as a bird :) 

Then you set_conn for the right element, set_conn ..( qr, 5632', ...) 

If you have no edges created in advance, you are done. 

If edge 25 was created in advance, you may want to also create edge 2'5 (with create_element method, or with get_adjacent methods, create flag true). 

Node 2' will not be connected at all with element on the left (4521), and node 2 will not be connected to element on the right. You will have a crack! 

Another note: in the above example, let's say that Range elems has 2 elements, those initial 2 quads 
when you request edges adjacent to them with a code like 
result = moab->get_adjacencies(elems,1,true,elems_edges,MBInterface::UNION); 
elem_edges will have exactly 7 edges: 14, 12, 45, 25, 56, 36, 23 
So edge 25 is not repeated, there is only one edge that connects 2 and 5. But it will be adjacent to both elements. Maybe your misunderstanding starts from here. 

At the same time, if you call that code after setting new connectivity on the element on the right, your elem_edges will have exactly 8 edges! 
(because you will have an extra edge, 2'5) 


You may want to look over "merge" methods. It is, in a way, the opposite of what you try to achieve, but it may give you a better understanding. 


I hope this helps. 

Iulian 


-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/moab-dev/attachments/20120504/6fd7b69f/attachment.htm>