nschloe / pygmsh

:spider_web: Gmsh for Python
GNU General Public License v3.0
843 stars 161 forks source link

Getting structured square cells in mesh #155

Closed tady57 closed 6 years ago

tady57 commented 6 years ago

Hi,

this is my python code (the important part): `with open('proba.txt', 'r') as myfile: for line in myfile: if 'Surface' in line:

        Surface.append(find_id(line))
        geom.add_raw_code(line)

    else:
        geom.add_raw_code(line)

myfile.close()

geom.add_raw_code('Mesh.Algorithm=8;')

for i in Surface:

geom.add_raw_code('Recombine Surface {%s}=0;' %i)

points, cells, point_data, cell_data, field_data = pygmsh.generate_mesh(geom)

meshio.write('test8.vtu', points, cells, cell_data=cell_data)

` And the pictures of mesh I get: nova4

nova5

The problem is that I have to have all squares or rectangles in my mesh and I am not able to get squares and rectangles at some surfaces(picture above). Do you have any advice how to repare that, which (py)gmsh command to use?

nschloe commented 6 years ago

Once again, let me refer you to https://guides.github.com/features/mastering-markdown/ for how to style code in GitHub markdown.

nschloe commented 6 years ago

Also, please read https://stackoverflow.com/help/mcve on how to create a verifiable example. Feel free to reopen with a code example that can be verified without proba.txt.

tady57 commented 6 years ago

proba.txt Here is the proba.txt file. It is basically copy paste .geo file with only points, lines and surfaces in it (which I get from gmsh). This all is working on small examples, but with huge geometry it does not. Here is the whole python code:

import pygmsh
import numpy as np
import meshio
from math import pow
import string
import sys
import traceback
import helpers
import scipy

def find_id(s):
    string_helper=''
    if'Surface' in s:
        s=s[8:]
        for i in s:
            if not i==')':
                string_helper=string_helper+i
            else:break
    string_helper=int(string_helper)        
    return string_helper

geom = pygmsh.built_in.Geometry()

Surface=[]

with open('proba.txt', 'r') as myfile:
    for line in myfile:
        if 'Surface' in line:

            Surface.append(find_id(line))
            geom.add_raw_code(line)

        else:
            geom.add_raw_code(line)

myfile.close() 

geom.add_raw_code('Mesh.Algorithm=8;')

for i in Surface:

    geom.add_raw_code('Recombine Surface {%s}=0;' %i)

points, cells, point_data, cell_data, field_data = pygmsh.generate_mesh(geom)

meshio.write('test8.vtu', points, cells, cell_data=cell_data)
nschloe commented 6 years ago

[...] without proba.txt.

Here is the proba.txt file.

You're attaching a 2000-line data file? Good luck expecting anyone to dig through that. I certainly won't.

tady57 commented 6 years ago

haha of course not. Nobody expecting you to go through that. It is fixed file with data from my professor. It must not be changed and it is here just to put geometry inside python code through geom.add_raw_code function. I put it here so you or someone else could compile my code and see the mesh. The problem is that on small examples it is very easy to get structed square mesh so my code's logic work perfectly. Maybe I should asked question like this: With fixed geometry and huge data what can I do to fix mesh - which algorithms to use, how to seperate surfaces so one algorithm works at one surface and other algorithm at some other surface? I have the problem with the exact surface (from picture above). It creates some funny deformed triangular / quadrilateral mesh. Thank you for you time, I am sorry for my imprecision.

nschloe commented 6 years ago

The problem is that on small examples it is very easy to get structed square mesh so my code's logic work perfectly.

Well, if it works for small meshes but doesn't for "huge", then there must be one smallest mesh for which is it doesn't. It's your job then to isolate this case. And frankly speaking: I don't think it depends of the size of the mesh.

nschloe commented 6 years ago

Also, I got to say that you're using pygmsh very weirdly. I mean: You already have a geo-file, right? Why don't you just create the mesh with it? Instead of your whole code above, you just append a bunch of Recombine Surfaces your geo file and be done with it.

PhilipHB commented 6 years ago

Maybe gmsh -check geometry.geo helps you.

tady57 commented 6 years ago

@nschloe Yes I know it is very easy way, but my professor gave me this .geo file and insists I have to use python

@phbphy Thank you, I will try that

nschloe commented 6 years ago

Yes I know it is very easy way, but my professor gave me this .geo file and insists I have to use python

I would be very interested to know why he makes you do that. To me, your code examples make no sense in pygmsh whatsoever.

tady57 commented 6 years ago

Ok, I will talk to him and explane situation. I will let you know the answer.

Thank you and have a nice evening

tady57 commented 6 years ago

Hi, I talked with him and he was not very interested in how to use pygmsh properly-the only condition is to use pygmsh. I am new in all this so its hard to just find solution. But to have meaningfull code I decided to create my own geometry with data given to me. On this example I am not able to create quad mesh on the cube surface. How can I change geom.add_raw_code('Recombine Surface {%s};' %cube) to work?

from helpers import compute_volume
import pygmsh
import numpy as np
import meshio
from math import pow
import string
import sys
import traceback
import helpers
import scipy

def test():
    lcar=pow(10,22)

    geom = pygmsh.built_in.Geometry()
    cube=geom.add_box(0, 1200, 0, 1200, 0, 1000, lcar)
    geom.add_raw_code('Recombine Surface {%s};' %cube)

    ref = 1.0
    points, cells, _, _, _ = pygmsh.generate_mesh(geom)
    #assert abs(compute_volume(points, cells) - ref) < 1.0e-2 * ref
    return points, cells

if __name__ == '__main__':
    import meshio
meshio.write('cube.vtu', *test())
nschloe commented 6 years ago

How can I change geom.add_raw_code('Recombine Surface {%s};' %cube) to work?

This problem has nothing to do with pygmsh, more with gmsh itself. Ask them, perhaps they'll be able to help out. https://gitlab.onelab.info/gmsh/gmsh

I talked with him and he was not very interested in how to use pygmsh properly-the only condition is to use pygmsh.

Wait. The head developer of the software tells you that this isn't how to do it, and he still insists on using it this way? Please give me his name and his email address. (You can send it to nico.schloemer@gmail.com.)

tady57 commented 6 years ago

Ok I will ask gmsh.

Thank you