test1.di

来自「算断裂的」· DI 代码 · 共 185 行

Thu Aug 26 14:52:32 EDT 1999
-----> ## QMG test 1


Thu Aug 26 14:52:32 EDT 1999
-----> ## This is an approximation to the mesh on the cover of Claes 


Thu Aug 26 14:52:32 EDT 1999
-----> ## Johnson's book ``Numerical solutions of partial differential


Thu Aug 26 14:52:32 EDT 1999
-----> ## equations by the finite element method'' , Cambridge U. Press 1987.


Thu Aug 26 14:52:32 EDT 1999
-----> ## It has greater refinement near the origin.


Thu Aug 26 14:52:32 EDT 1999
-----> 
global PI


Thu Aug 26 14:52:32 EDT 1999
-----> global interactive


Thu Aug 26 14:52:32 EDT 1999
-----> set wedgeangle [expr $PI/6.0]

0.52359877559829837
Thu Aug 26 14:52:32 EDT 1999
-----> 

## Make a circular approximation with 4 Bezier curves from pi/6 to 2pi.


Thu Aug 26 14:52:32 EDT 1999
-----> 
gmset {verts codes} [gm_circ_approx $wedgeangle [expr 2 * $PI] 4]


Thu Aug 26 14:52:32 EDT 1999
-----> 
## Put one more vertex to start the wedge.


Thu Aug 26 14:52:32 EDT 1999
-----> 
gmset newverts [concat {{0 0}} $verts]


Thu Aug 26 14:52:32 EDT 1999
-----> gmset newcodes [concat 0 $codes]


Thu Aug 26 14:52:32 EDT 1999
-----> gmset domain [gm_cpoly $newverts $newcodes]


Thu Aug 26 14:52:32 EDT 1999
-----> 

## The size function is of the form h0*r^p+f, where p is determined


Thu Aug 26 14:52:32 EDT 1999
-----> ## by (4.31) in Johnson, and f is determined so that the H1 error on


Thu Aug 26 14:52:32 EDT 1999
-----> ## the innermost triangles matches the H1 error everywhere else.


Thu Aug 26 14:52:32 EDT 1999
-----> 
set h0 0.2

0.2
Thu Aug 26 14:52:32 EDT 1999
-----> set gamma [expr $PI / ( 2*$PI - $wedgeangle)]

0.54545454545454541
Thu Aug 26 14:52:32 EDT 1999
-----> set f [expr pow($h0,(1.0/$gamma))]

0.052306419440473213
Thu Aug 26 14:52:32 EDT 1999
-----> set p [expr 1.0 - $gamma]

0.45454545454545459
Thu Aug 26 14:52:32 EDT 1999
-----> 
set sizectl "(formula ($h0 * pow(%0*%0+%1*%1, $p/2.0) + $f))"

(formula (0.2 * pow(%0*%0+%1*%1, 0.45454545454545459/2.0) + 0.052306419440473213))
Thu Aug 26 14:52:32 EDT 1999
-----> 
set show 0

0
Thu Aug 26 14:52:32 EDT 1999
-----> if {[llength [info globals interactive]]} {
    puts "size control string = $sizectl"
    set show 1
}


Thu Aug 26 14:52:32 EDT 1999
-----> 

gmset mesh [gmmeshgen $domain size $sizectl show $show]


Thu Aug 26 14:52:39 EDT 1999
-----> gmset {numvtx numelt} [gmmeshsize $mesh]


Thu Aug 26 14:52:39 EDT 1999
-----> set asp [gmchecktri $domain $mesh]

Maximum aspect ratio =        10.2397 achieved in
simplex #359 of topological entity mregion (2:0) which has vertices 13 11 230 
Maximum global side length =  0.200806
Minimum global altitude =     0.00727931
Number of nodes = 259 number of elements = 447
10.239729237340125
Thu Aug 26 14:52:40 EDT 1999
-----> 
if {[llength [info globals aspprod]]} {
    global aspprod
    set aspprod [expr $aspprod * $asp]
    global meshsizesum
    set meshsizesum [expr $meshsizesum + $numvtx]
}

259
Thu Aug 26 14:52:40 EDT 1999
-----> 
if {[llength [info globals interactive]]} {
    gmviz $mesh white
}


Thu Aug 26 14:52:40 EDT 1999
-----> 
# ------------------------------------------------------------------


Thu Aug 26 14:52:40 EDT 1999
-----> # Copyright (c) 1999 by Cornell University.  All rights reserved


Thu Aug 26 14:52:40 EDT 1999
-----> # See the accompanying file 'Copyright' for authorship information,


Thu Aug 26 14:52:40 EDT 1999
-----> # the terms of the license governing this software, and disclaimers


Thu Aug 26 14:52:40 EDT 1999
-----> # concerning this software.


Thu Aug 26 14:52:40 EDT 1999
-----> # ------------------------------------------------------------------


Thu Aug 26 14:52:40 EDT 1999
-----> # This file is part of the QMG software.  


Thu Aug 26 14:52:40 EDT 1999
-----> # Version 2.0 of QMG, release date RELDATE.


Thu Aug 26 14:52:40 EDT 1999
-----> # ------------------------------------------------------------------


Thu Aug 26 14:52:40 EDT 1999