[CIG-SEISMO] Negative Jacobian in specfem2d
Jami Johnson
jami.johnson at auckland.ac.nz
Wed Aug 3 18:58:10 PDT 2016
Hi Dimitri and Paul,
Thank you for your replies! I also have run the mesh and received no
errors when I run xcheck_quality_external_mesh. After another test, I
found that when I run SPECFEM2D with a mesh size of 0.00013, it runs
without error, but when I decrease the element size to 0.00003, I get the
negative Jacobian error (even when xcheck_quality_external_mesh shows no
error). Is there a resolution limit? In my parameter file, I use
deltat=5d-10 and velocity of 1500 m/s.
Thanks again for your help.
Cheers,
Jami
On Thu, Aug 4, 2016 at 3:14 AM cristini <cristini at lma.cnrs-mrs.fr> wrote:
> Hi Jami,
>
> I did a quick test and did not find any problem. This is probably a Gmsh
> issue that you should fix easily by upgrading to a more recent version.
>
> I took you .geo file and got a very good mesh without any negative
> jacobian.
>
> Thanks,
>
> Best wishes
>
> Paul
>
> PS: the histogram I got
>
> histogram of skewness (0. good - 1. bad):
>
> 0.0000000E+00 - 5.0000001E-02 984857 95.14644 %
> 5.0000001E-02 - 0.1000000 40000 3.864376 %
> 0.1000000 - 0.1500000 2983 0.2881858 %
> 0.1500000 - 0.2000000 2547 0.2460641 %
> 0.2000000 - 0.2500000 1288 0.1244329 %
> 0.2500000 - 0.3000000 1358 0.1311956 %
> 0.3000000 - 0.3500000 597 5.7675809E-02 %
> 0.3500000 - 0.4000000 595 5.7482593E-02 %
> 0.4000000 - 0.4500000 619 5.9801217E-02 %
> 0.4500000 - 0.5000000 189 1.8259175E-02 %
> 0.5000000 - 0.5500000 61 5.8931732E-03 %
> 0.5500000 - 0.6000000 1 9.6609394E-05 %
> 0.6000000 - 0.6500000 0 0.0000000E+00 %
> 0.6500000 - 0.7000000 1 9.6609394E-05 %
> 0.7000000 - 0.7500000 0 0.0000000E+00 %
> 0.7500000 - 0.8000000 0 0.0000000E+00 %
> 0.8000000 - 0.8500000 0 0.0000000E+00 %
> 0.8500000 - 0.9000000 0 0.0000000E+00 %
> 0.9000000 - 0.9500000 0 0.0000000E+00 %
> 0.9500000 - 1.000000 0 0.0000000E+00 %
>
> total number of elements = 1035096
>
>
>
>
> On 03/08/2016 15:43, Dimitri Komatitsch wrote:
> >
> > Hi Jami,
> >
> > Thanks! Probably mesh quality that is not good enough, two elements
> > that are almost flat. Try to see the quality histogram in Gmsh to see
> > if these two elements are poor (in which case even if you manage to
> > run that mesh in SPECFEM2D it will likely blow up). There is also a
> > tool for that in SPECFEM2D:
> >
> > src/auxiliaries/check_quality_external_mesh.f90
> >
> > If it turns out to be a bug in the code (unlikely, but worth checking)
> > please let us know.
> >
> > Thanks,
> > Best wishes,
> >
> > Dimitri.
> >
> > On 03/08/2016 10:01, Jami Johnson wrote:
> >> Hello,
> >> I have created a very dense mesh with Gmsh, and when I check for the
> >> values of the Jacobian in Gmsh they are all positive. However, when I
> >> run it in specfem2d, an error is thrown due to a negative Jacobian. It
> >> seems to be a very small percentage of my elements causing a problem. I
> >> believe all of my line loops are oriented correctly. Are there any
> >> other potential causes of a negative Jacobian?
> >>
> >> Error in specfem2d:
> >> --------------------------
> >> 2 elements have a negative Jacobian, out of 1313871
> >>
> >> geo file:
> >> -----------
> >> lc=0.00003;
> >> xmax=0.015;
> >> xmin=0.0;
> >> ymin=-0.015;
> >> ymax=0.0;
> >> Point(1) = {xmax,ymax , 0, lc};
> >> Point(2) = {xmax, ymin, 0, lc};
> >> Point(3) = {xmin, ymin, 0, lc};
> >> Point(4) = {xmin, ymax, 0, lc};
> >> Line(1) = {1, 4};
> >> Line(2) = {4, 3};
> >> Line(3) = {3, 2};
> >> Line(4) = {2, 1};
> >> R=0.00025;
> >> Centx=0.0075;
> >> Centy=-0.005;
> >> Point(5) = {Centx, R+Centy, 0, lc};
> >> Point(6) = {R+Centx, Centy, 0, lc};
> >> Point(7) = {-R+Centx, Centy, 0, lc};
> >> Point(8) = {Centx, -R+Centy, 0, lc};
> >> Point(9) = {Centx, Centy, 0, lc};
> >> Circle(5) = {5, 9, 7};
> >> Circle(6) = {7, 9, 8};
> >> Circle(7) = {8, 9, 6};
> >> Circle(8) = {6, 9, 5};
> >> Line Loop(9) = {5, 6, 7, 8};
> >> Line Loop(10) = {1, 2, 3, 4};
> >>
> >> Plane Surface(10) = {10, 9};
> >> Plane Surface(11) = {9};
> >>
> >> Recombine Surface{10,11};
> >> Mesh.SubdivisionAlgorithm=1;
> >> Physical Line("Top") = {1};
> >> Physical Line("Left") = {2};
> >> Physical Line("Bottom") = {3};
> >> Physical Line("Right") = {4};
> >> Physical Surface("M1") = {10};
> >> Physical Surface("M2") = {11};
> >>
> >>
> >>
> >>
> >> _______________________________________________
> >> CIG-SEISMO mailing list
> >> CIG-SEISMO at geodynamics.org
> >> http://lists.geodynamics.org/cgi-bin/mailman/listinfo/cig-seismo
> >>
> >
>
> --
>
> *--------------------------------------------------------------------------*
> Warning !!! New address !! Attention !!! nouvelle adresse
> !!
>
> *--------------------------------------------------------------------------*
> Paul CRISTINI Charge de Recherche CNRS (HDR) / CNRS Research Scientist
>
> *--------------------------------------------------------------------------*
> CNRS - Laboratoire de Mecanique et d'Acoustique (UPR 7051)
> Bureau 119
> 4 Impasse Nikola Tesla, CS 40006, F-13453 Marseille Cedex 13 - France
> mailto: cristini at lma.cnrs-mrs.fr
> phone number: +33 (0)4 84 52 42 51
> http://www.lma.cnrs-mrs.fr/
>
> *--------------------------------------------------------------------------*
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.geodynamics.org/pipermail/cig-seismo/attachments/20160804/6ce7d96c/attachment.html>
More information about the CIG-SEISMO
mailing list