-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathintpts.f90
668 lines (658 loc) · 40.1 KB
/
intpts.f90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
subroutine Gauss1D(p,x,w)
implicit none
integer, intent(in) :: p ! p is integration order! not solution order
real(8), intent(out), dimension((p+2)/2) :: x,w ! integer division to apply floor, since it's really floor((p+1)/2)
if (p <= 1) then
! Order 1d0 Gauss-Legendre points
x = (/ &
0.500000000000000d0 &
/)
w = (/ &
1.000000000000000d0 &
/)
elseif (p <= 3) then
! Order 3d0 Gauss-Legendre points
x = (/ &
0.211324865405187d0, 0.788675134594813d0 &
/)
w = (/ &
0.500000000000000d0, 0.500000000000000d0 &
/)
elseif (p <= 5) then
! Order 5d0 Gauss-Legendre points
x = (/ &
0.112701665379258d0, 0.500000000000000d0, 0.887298334620742d0 &
/)
w = (/ &
0.277777777777778d0, 0.444444444444444d0, 0.277777777777778d0 &
/)
elseif (p <= 7) then
! Order 7d0 Gauss-Legendre points
x = (/ &
0.069431844202974d0, 0.330009478207572d0, 0.669990521792428d0, 0.930568155797026d0 &
/)
w = (/ &
0.173927422568727d0, 0.326072577431273d0, 0.326072577431273d0, 0.173927422568727d0 &
/)
elseif (p <= 9) then
! Order 9d0 Gauss-Legendre points
x = (/ &
0.046910077030668d0, 0.230765344947158d0, 0.500000000000000d0, 0.769234655052841d0, &
0.953089922969332d0 &
/)
w = (/ &
0.118463442528095d0, 0.239314335249683d0, 0.284444444444444d0, 0.239314335249683d0, &
0.118463442528095d0 &
/)
elseif (p <= 11) then
! Order 11d0 Gauss-Legendre points
x = (/ &
0.033765242898424d0, 0.169395306766868d0, 0.380690406958402d0, 0.619309593041598d0, &
0.830604693233132d0, 0.966234757101576d0 &
/)
w = (/ &
0.085662246189585d0, 0.180380786524069d0, 0.233956967286345d0, 0.233956967286345d0, &
0.180380786524069d0, 0.085662246189585d0 &
/)
elseif (p <= 13) then
! Order 13d0 Gauss-Legendre points
x = (/ &
0.025446043828621d0, 0.129234407200303d0, 0.297077424311301d0, 0.500000000000000d0, &
0.702922575688699d0, 0.870765592799697d0, 0.974553956171379d0 &
/)
w = (/ &
0.064742483084435d0, 0.139852695744638d0, 0.190915025252560d0, 0.208979591836735d0, &
0.190915025252560d0, 0.139852695744638d0, 0.064742483084435d0 &
/)
elseif (p <= 15) then
! Order 15d0 Gauss-Legendre points
x = (/ &
0.019855071751232d0, 0.101666761293187d0, 0.237233795041836d0, 0.408282678752175d0, &
0.591717321247825d0, 0.762766204958164d0, 0.898333238706813d0, 0.980144928248768d0 &
/)
w = (/ &
0.050614268145188d0, 0.111190517226687d0, 0.156853322938944d0, 0.181341891689181d0, &
0.181341891689181d0, 0.156853322938944d0, 0.111190517226687d0, 0.050614268145188d0 &
/)
elseif (p <= 17) then
! Order 17d0 Gauss-Legendre points
x = (/ &
0.015919880246187d0, 0.081984446336682d0, 0.193314283649705d0, 0.337873288298096d0, &
0.500000000000000d0, 0.662126711701905d0, 0.806685716350295d0, 0.918015553663318d0, &
0.984080119753813d0 &
/)
w = (/ &
0.040637194180787d0, 0.090324080347429d0, 0.130305348201468d0, 0.156173538520001d0, &
0.165119677500630d0, 0.156173538520001d0, 0.130305348201468d0, 0.090324080347429d0, &
0.040637194180787d0 &
/)
elseif (p <= 19) then
! Order 19d0 Gauss-Legendre points
x = (/ &
0.013046735741414d0, 0.067468316655508d0, 0.160295215850488d0, 0.283302302935376d0, &
0.425562830509184d0, 0.574437169490816d0, 0.716697697064624d0, 0.839704784149512d0, &
0.932531683344492d0, 0.986953264258586d0 &
/)
w = (/ &
0.033335672154344d0, 0.074725674575290d0, 0.109543181257991d0, 0.134633359654998d0, &
0.147762112357376d0, 0.147762112357376d0, 0.134633359654998d0, 0.109543181257991d0, &
0.074725674575290d0, 0.033335672154344d0 &
/)
elseif (p <= 21) then
! Order 21d0 Gauss-Legendre points
x = (/ &
0.010885670926972d0, 0.056468700115952d0, 0.134923997212975d0, 0.240451935396594d0, &
0.365228422023828d0, 0.500000000000000d0, 0.634771577976172d0, 0.759548064603406d0, &
0.865076002787025d0, 0.943531299884048d0, 0.989114329073028d0 &
/)
w = (/ &
0.027834283558087d0, 0.062790184732452d0, 0.093145105463867d0, 0.116596882295995d0, &
0.131402272255123d0, 0.136462543388950d0, 0.131402272255123d0, 0.116596882295995d0, &
0.093145105463867d0, 0.062790184732452d0, 0.027834283558087d0 &
/)
elseif (p <= 23) then
! Order 23d0 Gauss-Legendre points
x = (/ &
0.009219682876640d0, 0.047941371814763d0, 0.115048662902848d0, 0.206341022856691d0, &
0.316084250500910d0, 0.437383295744266d0, 0.562616704255734d0, 0.683915749499090d0, &
0.793658977143309d0, 0.884951337097152d0, 0.952058628185237d0, 0.990780317123360d0 &
/)
w = (/ &
0.023587668193256d0, 0.053469662997659d0, 0.080039164271673d0, 0.101583713361533d0, &
0.116746268269177d0, 0.124573522906701d0, 0.124573522906701d0, 0.116746268269177d0, &
0.101583713361533d0, 0.080039164271673d0, 0.053469662997659d0, 0.023587668193256d0 &
/)
elseif (p <= 25) then
! Order 25d0 Gauss-Legendre points
x = (/ &
0.007908472640706d0, 0.041200800388511d0, 0.099210954633345d0, 0.178825330279830d0, &
0.275753624481777d0, 0.384770842022433d0, 0.500000000000000d0, 0.615229157977567d0, &
0.724246375518223d0, 0.821174669720170d0, 0.900789045366655d0, 0.958799199611489d0, &
0.992091527359294d0 &
/)
w = (/ &
0.020242002382658d0, 0.046060749918864d0, 0.069436755109894d0, 0.089072990380973d0, &
0.103908023768444d0, 0.113141590131449d0, 0.116275776615437d0, 0.113141590131449d0, &
0.103908023768444d0, 0.089072990380973d0, 0.069436755109894d0, 0.046060749918864d0, &
0.020242002382658d0 &
/)
elseif (p <= 27) then
! Order 27d0 Gauss-Legendre points
x = (/ &
0.006858095651594d0, 0.035782558168213d0, 0.086399342465117d0, 0.156353547594157d0, &
0.242375681820923d0, 0.340443815536055d0, 0.445972525646328d0, 0.554027474353672d0, &
0.659556184463945d0, 0.757624318179077d0, 0.843646452405843d0, 0.913600657534883d0, &
0.964217441831787d0, 0.993141904348406d0 &
/)
w = (/ &
0.017559730165876d0, 0.040079043579880d0, 0.060759285343952d0, 0.078601583579097d0, &
0.092769198738969d0, 0.102599231860648d0, 0.107631926731579d0, 0.107631926731579d0, &
0.102599231860648d0, 0.092769198738969d0, 0.078601583579097d0, 0.060759285343952d0, &
0.040079043579880d0, 0.017559730165876d0 &
/)
elseif (p <= 29) then
! Order 29d0 Gauss-Legendre points
x = (/ &
0.006003740989757d0, 0.031363303799647d0, 0.075896708294786d0, 0.137791134319915d0, &
0.214513913695731d0, 0.302924326461218d0, 0.399402953001283d0, 0.500000000000000d0, &
0.600597046998717d0, 0.697075673538782d0, 0.785486086304269d0, 0.862208865680085d0, &
0.924103291705214d0, 0.968636696200353d0, 0.993996259010243d0 &
/)
w = (/ &
0.015376620998059d0, 0.035183023744054d0, 0.053579610233586d0, 0.069785338963077d0, &
0.083134602908497d0, 0.093080500007781d0, 0.099215742663556d0, 0.101289120962781d0, &
0.099215742663556d0, 0.093080500007781d0, 0.083134602908497d0, 0.069785338963077d0, &
0.053579610233586d0, 0.035183023744054d0, 0.015376620998059d0 &
/)
elseif (p <= 31) then
! Order 31d0 Gauss-Legendre points
x = (/ &
0.005299532504175d0, 0.027712488463384d0, 0.067184398806084d0, 0.122297795822498d0, &
0.191061877798678d0, 0.270991611171386d0, 0.359198224610371d0, 0.452493745081181d0, &
0.547506254918819d0, 0.640801775389629d0, 0.729008388828614d0, 0.808938122201322d0, &
0.877702204177502d0, 0.932815601193916d0, 0.972287511536616d0, 0.994700467495825d0 &
/)
w = (/ &
0.013576229705877d0, 0.031126761969324d0, 0.047579255841246d0, 0.062314485627767d0, &
0.074797994408288d0, 0.084578259697501d0, 0.091301707522462d0, 0.094725305227534d0, &
0.094725305227534d0, 0.091301707522462d0, 0.084578259697501d0, 0.074797994408288d0, &
0.062314485627767d0, 0.047579255841246d0, 0.031126761969324d0, 0.013576229705877d0 &
/)
elseif (p <= 33) then
! Order 33d0 Gauss-Legendre points
x = (/ &
0.004712262342791d0, 0.024662239115616d0, 0.059880423136507d0, 0.109242998051599d0, &
0.171164420391655d0, 0.243654731456762d0, 0.324384118273062d0, 0.410757909252076d0, &
0.500000000000000d0, 0.589242090747924d0, 0.675615881726938d0, 0.756345268543238d0, &
0.828835579608345d0, 0.890757001948401d0, 0.940119576863493d0, 0.975337760884384d0, &
0.995287737657209d0 &
/)
w = (/ &
0.012074151434274d0, 0.027729764686994d0, 0.042518074158590d0, 0.055941923596702d0, &
0.067568184234263d0, 0.077022880538405d0, 0.084002051078225d0, 0.088281352683496d0, &
0.089723235178103d0, 0.088281352683496d0, 0.084002051078225d0, 0.077022880538405d0, &
0.067568184234263d0, 0.055941923596702d0, 0.042518074158590d0, 0.027729764686994d0, &
0.012074151434274d0 &
/)
elseif (p <= 35) then
! Order 35d0 Gauss-Legendre points
x = (/ &
0.004217415789535d0, 0.022088025214301d0, 0.053698766751222d0, 0.098147520513738d0, &
0.154156478469823d0, 0.220114584463026d0, 0.294124419268579d0, 0.374056887154247d0, &
0.457612493479132d0, 0.542387506520868d0, 0.625943112845753d0, 0.705875580731421d0, &
0.779885415536974d0, 0.845843521530177d0, 0.901852479486262d0, 0.946301233248778d0, &
0.977911974785699d0, 0.995782584210466d0 &
/)
w = (/ &
0.010808006763242d0, 0.024857274447485d0, 0.038212865127445d0, 0.050471022053144d0, &
0.061277603355739d0, 0.070321457335325d0, 0.077342337563133d0, 0.082138241872916d0, &
0.084571191481572d0, 0.084571191481572d0, 0.082138241872916d0, 0.077342337563133d0, &
0.070321457335325d0, 0.061277603355739d0, 0.050471022053144d0, 0.038212865127445d0, &
0.024857274447485d0, 0.010808006763242d0 &
/)
elseif (p <= 37) then
! Order 37d0 Gauss-Legendre points
x = (/ &
0.003796578078208d0, 0.019895923932585d0, 0.048422048192591d0, 0.088642671731429d0, &
0.139516911332385d0, 0.199727347669160d0, 0.267714629312020d0, 0.341717950018185d0, &
0.419820677179887d0, 0.500000000000000d0, 0.580179322820113d0, 0.658282049981815d0, &
0.732285370687980d0, 0.800272652330841d0, 0.860483088667615d0, 0.911357328268571d0, &
0.951577951807409d0, 0.980104076067415d0, 0.996203421921792d0 &
/)
w = (/ &
0.009730894114863d0, 0.022407113382850d0, 0.034522271368821d0, 0.045745010811225d0, &
0.055783322773667d0, 0.064376981269668d0, 0.071303351086803d0, 0.076383021032930d0, &
0.079484421696977d0, 0.080527224924392d0, 0.079484421696977d0, 0.076383021032930d0, &
0.071303351086803d0, 0.064376981269668d0, 0.055783322773667d0, 0.045745010811225d0, &
0.034522271368821d0, 0.022407113382850d0, 0.009730894114863d0 &
/)
elseif (p <= 39) then
! Order 39d0 Gauss-Legendre points
x = (/ &
0.003435700407453d0, 0.018014036361043d0, 0.043882785874337d0, 0.080441514088891d0, &
0.126834046769925d0, 0.181973159636742d0, 0.244566499024586d0, 0.313146955642290d0, &
0.386107074429177d0, 0.461736739433251d0, 0.538263260566749d0, 0.613892925570823d0, &
0.686853044357710d0, 0.755433500975414d0, 0.818026840363258d0, 0.873165953230075d0, &
0.919558485911109d0, 0.956117214125663d0, 0.981985963638957d0, 0.996564299592547d0 &
/)
w = (/ &
0.008807003569576d0, 0.020300714900193d0, 0.031336024167055d0, 0.041638370788352d0, &
0.050965059908620d0, 0.059097265980759d0, 0.065844319224588d0, 0.071048054659191d0, &
0.074586493236302d0, 0.076376693565363d0, 0.076376693565363d0, 0.074586493236302d0, &
0.071048054659191d0, 0.065844319224588d0, 0.059097265980759d0, 0.050965059908620d0, &
0.041638370788352d0, 0.031336024167055d0, 0.020300714900193d0, 0.008807003569576d0 &
/)
endif
end subroutine Gauss1D
subroutine Gauss2D_pre(p,N)
implicit none
integer, intent(in) :: p ! again p is integration order, i.e. g
integer, intent(out) :: N
if (p <= 1) then
N = 1
elseif (p == 2) then
N = 3
elseif (p == 3) then
N = 4
elseif (p == 4) then
N = 6
elseif (p == 5) then
N = 7
elseif (p == 6) then
N = 12
elseif (p == 7) then
N = 13
elseif (p == 8) then
N = 16
elseif (p == 9) then
N = 19
elseif (p == 10) then
N = 25
elseif (p <= 12) then
N = 33
elseif (p == 13) then
N = 37
elseif (p == 14) then
N = 42
elseif (p <= 17) then
N = 61
elseif (p <= 19) then
N = 73
endif
end subroutine Gauss2D_pre
subroutine Gauss2D(p,N,xy,w)
implicit none
integer, intent(in) :: p
integer, intent(in) :: N
real(8), intent(out), dimension(N,2) :: xy
real(8), intent(out), dimension(N) :: w
real(8), dimension(2*N) :: x
if (p <= 1) then
! Order 1d0 Dunavant Points
x = (/ &
0.333333333333333d0, 0.333333333333333d0 &
/)
w = (/ &
0.500000000000000d0 &
/)
elseif (p == 2) then
! Order 2d0 Dunavant Points
x = (/ &
0.666666666666667d0, 0.166666666666667d0, 0.166666666666667d0, 0.166666666666667d0, &
0.166666666666667d0, 0.666666666666667d0 &
/)
w = (/ &
0.166666666666666d0, 0.166666666666666d0, 0.166666666666666d0 &
/)
elseif (p == 3) then
! Order 3d0 Dunavant Points
x = (/ &
0.333333333333333d0, 0.333333333333333d0, 0.600000000000000d0, 0.200000000000000d0, &
0.200000000000000d0, 0.200000000000000d0, 0.200000000000000d0, 0.600000000000000d0 &
/)
w = (/ &
-0.281250000000000d0, 0.260416666666667d0, 0.260416666666667d0, 0.260416666666667d0 &
/)
elseif (p == 4) then
! Order 4d0 Dunavant Points
x = (/ &
0.108103018168070d0, 0.445948490915965d0, 0.445948490915965d0, 0.445948490915965d0, &
0.445948490915965d0, 0.108103018168070d0, 0.816847572980459d0, 0.091576213509771d0, &
0.091576213509771d0, 0.091576213509771d0, 0.091576213509771d0, 0.816847572980459d0 &
/)
w = (/ &
0.111690794839005d0, 0.111690794839005d0, 0.111690794839005d0, 0.054975871827661d0, &
0.054975871827661d0, 0.054975871827661d0 &
/)
elseif (p == 5) then
! Order 5d0 Dunavant Points
x = (/ &
0.333333333333333d0, 0.333333333333333d0, 0.059715871789770d0, 0.470142064105115d0, &
0.470142064105115d0, 0.470142064105115d0, 0.470142064105115d0, 0.059715871789770d0, &
0.797426985353087d0, 0.101286507323456d0, 0.101286507323456d0, 0.101286507323456d0, &
0.101286507323456d0, 0.797426985353087d0 &
/)
w = (/ &
0.112500000000000d0, 0.066197076394253d0, 0.066197076394253d0, 0.066197076394253d0, &
0.062969590272414d0, 0.062969590272414d0, 0.062969590272414d0 &
/)
elseif (p == 6) then
! Order 6d0 Dunavant Points
x = (/ &
0.501426509658179d0, 0.249286745170910d0, 0.249286745170910d0, 0.249286745170910d0, &
0.249286745170910d0, 0.501426509658179d0, 0.873821971016996d0, 0.063089014491502d0, &
0.063089014491502d0, 0.063089014491502d0, 0.063089014491502d0, 0.873821971016996d0, &
0.053145049844817d0, 0.310352451033784d0, 0.310352451033784d0, 0.636502499121399d0, &
0.636502499121399d0, 0.053145049844817d0, 0.310352451033784d0, 0.053145049844817d0, &
0.636502499121399d0, 0.310352451033784d0, 0.053145049844817d0, 0.636502499121399d0 &
/)
w = (/ &
0.058393137863189d0, 0.058393137863189d0, 0.058393137863189d0, 0.025422453185103d0, &
0.025422453185103d0, 0.025422453185103d0, 0.041425537809187d0, 0.041425537809187d0, &
0.041425537809187d0, 0.041425537809187d0, 0.041425537809187d0, 0.041425537809187d0 &
/)
elseif (p == 7) then
! Order 7d0 Dunavant Points
x = (/ &
0.333333333333333d0, 0.333333333333333d0, 0.479308067841920d0, 0.260345966079040d0, &
0.260345966079040d0, 0.260345966079040d0, 0.260345966079040d0, 0.479308067841920d0, &
0.869739794195568d0, 0.065130102902216d0, 0.065130102902216d0, 0.065130102902216d0, &
0.065130102902216d0, 0.869739794195568d0, 0.048690315425316d0, 0.312865496004874d0, &
0.312865496004874d0, 0.638444188569810d0, 0.638444188569810d0, 0.048690315425316d0, &
0.312865496004874d0, 0.048690315425316d0, 0.638444188569810d0, 0.312865496004874d0, &
0.048690315425316d0, 0.638444188569810d0 &
/)
w = (/ &
-0.074785022233841d0, 0.087807628716604d0, 0.087807628716604d0, 0.087807628716604d0, &
0.026673617804419d0, 0.026673617804419d0, 0.026673617804419d0, 0.038556880445128d0, &
0.038556880445128d0, 0.038556880445128d0, 0.038556880445128d0, 0.038556880445128d0, &
0.038556880445128d0 &
/)
elseif (p == 8) then
! Order 8d0 Dunavant Points
x = (/ &
0.333333333333333d0, 0.333333333333333d0, 0.081414823414554d0, 0.459292588292723d0, &
0.459292588292723d0, 0.459292588292723d0, 0.459292588292723d0, 0.081414823414554d0, &
0.658861384496480d0, 0.170569307751760d0, 0.170569307751760d0, 0.170569307751760d0, &
0.170569307751760d0, 0.658861384496480d0, 0.898905543365938d0, 0.050547228317031d0, &
0.050547228317031d0, 0.050547228317031d0, 0.050547228317031d0, 0.898905543365938d0, &
0.008394777409958d0, 0.263112829634638d0, 0.263112829634638d0, 0.728492392955404d0, &
0.728492392955404d0, 0.008394777409958d0, 0.263112829634638d0, 0.008394777409958d0, &
0.728492392955404d0, 0.263112829634638d0, 0.008394777409958d0, 0.728492392955404d0 &
/)
w = (/ &
0.072157803838894d0, 0.047545817133642d0, 0.047545817133642d0, 0.047545817133642d0, &
0.051608685267359d0, 0.051608685267359d0, 0.051608685267359d0, 0.016229248811599d0, &
0.016229248811599d0, 0.016229248811599d0, 0.013615157087217d0, 0.013615157087217d0, &
0.013615157087217d0, 0.013615157087217d0, 0.013615157087217d0, 0.013615157087217d0 &
/)
elseif (p == 9) then
! Order 9d0 Dunavant Points
x = (/ &
0.333333333333333d0, 0.333333333333333d0, 0.020634961602525d0, 0.489682519198738d0, &
0.489682519198738d0, 0.489682519198738d0, 0.489682519198738d0, 0.020634961602525d0, &
0.125820817014127d0, 0.437089591492937d0, 0.437089591492937d0, 0.437089591492937d0, &
0.437089591492937d0, 0.125820817014127d0, 0.623592928761935d0, 0.188203535619033d0, &
0.188203535619033d0, 0.188203535619033d0, 0.188203535619033d0, 0.623592928761935d0, &
0.910540973211095d0, 0.044729513394453d0, 0.044729513394453d0, 0.044729513394453d0, &
0.044729513394453d0, 0.910540973211095d0, 0.036838412054736d0, 0.221962989160766d0, &
0.221962989160766d0, 0.741198598784498d0, 0.741198598784498d0, 0.036838412054736d0, &
0.221962989160766d0, 0.036838412054736d0, 0.741198598784498d0, 0.221962989160766d0, &
0.036838412054736d0, 0.741198598784498d0 &
/)
w = (/ &
0.048567898141400d0, 0.015667350113570d0, 0.015667350113570d0, 0.015667350113570d0, &
0.038913770502387d0, 0.038913770502387d0, 0.038913770502387d0, 0.039823869463605d0, &
0.039823869463605d0, 0.039823869463605d0, 0.012788837829349d0, 0.012788837829349d0, &
0.012788837829349d0, 0.021641769688645d0, 0.021641769688645d0, 0.021641769688645d0, &
0.021641769688645d0, 0.021641769688645d0, 0.021641769688645d0 &
/)
elseif (p == 10) then
! Order 10d0 Dunavant Points
x = (/ &
0.333333333333333d0, 0.333333333333333d0, 0.028844733232685d0, 0.485577633383657d0, &
0.485577633383657d0, 0.485577633383657d0, 0.485577633383657d0, 0.028844733232685d0, &
0.781036849029926d0, 0.109481575485037d0, 0.109481575485037d0, 0.109481575485037d0, &
0.109481575485037d0, 0.781036849029926d0, 0.141707219414880d0, 0.307939838764121d0, &
0.307939838764121d0, 0.550352941820999d0, 0.550352941820999d0, 0.141707219414880d0, &
0.307939838764121d0, 0.141707219414880d0, 0.550352941820999d0, 0.307939838764121d0, &
0.141707219414880d0, 0.550352941820999d0, 0.025003534762686d0, 0.246672560639903d0, &
0.246672560639903d0, 0.728323904597411d0, 0.728323904597411d0, 0.025003534762686d0, &
0.246672560639903d0, 0.025003534762686d0, 0.728323904597411d0, 0.246672560639903d0, &
0.025003534762686d0, 0.728323904597411d0, 0.009540815400299d0, 0.066803251012200d0, &
0.066803251012200d0, 0.923655933587500d0, 0.923655933587500d0, 0.009540815400299d0, &
0.066803251012200d0, 0.009540815400299d0, 0.923655933587500d0, 0.066803251012200d0, &
0.009540815400299d0, 0.923655933587500d0 &
/)
w = (/ &
0.045408995191377d0, 0.018362978878233d0, 0.018362978878233d0, 0.018362978878233d0, &
0.022660529717764d0, 0.022660529717764d0, 0.022660529717764d0, 0.036378958422710d0, &
0.036378958422710d0, 0.036378958422710d0, 0.036378958422710d0, 0.036378958422710d0, &
0.036378958422710d0, 0.014163621265528d0, 0.014163621265528d0, 0.014163621265528d0, &
0.014163621265528d0, 0.014163621265528d0, 0.014163621265528d0, 0.004710833481867d0, &
0.004710833481867d0, 0.004710833481867d0, 0.004710833481867d0, 0.004710833481867d0, &
0.004710833481867d0 &
/)
elseif (p <= 12) then
! Order 12d0 Dunavant Points
x = (/ &
0.023565220452390d0, 0.488217389773805d0, 0.488217389773805d0, 0.488217389773805d0, &
0.488217389773805d0, 0.023565220452390d0, 0.120551215411079d0, 0.439724392294460d0, &
0.439724392294460d0, 0.439724392294460d0, 0.439724392294460d0, 0.120551215411079d0, &
0.457579229975768d0, 0.271210385012116d0, 0.271210385012116d0, 0.271210385012116d0, &
0.271210385012116d0, 0.457579229975768d0, 0.744847708916828d0, 0.127576145541586d0, &
0.127576145541586d0, 0.127576145541586d0, 0.127576145541586d0, 0.744847708916828d0, &
0.957365299093579d0, 0.021317350453210d0, 0.021317350453210d0, 0.021317350453210d0, &
0.021317350453210d0, 0.957365299093579d0, 0.115343494534698d0, 0.275713269685514d0, &
0.275713269685514d0, 0.608943235779788d0, 0.608943235779788d0, 0.115343494534698d0, &
0.275713269685514d0, 0.115343494534698d0, 0.608943235779788d0, 0.275713269685514d0, &
0.115343494534698d0, 0.608943235779788d0, 0.022838332222257d0, 0.281325580989940d0, &
0.281325580989940d0, 0.695836086787803d0, 0.695836086787803d0, 0.022838332222257d0, &
0.281325580989940d0, 0.022838332222257d0, 0.695836086787803d0, 0.281325580989940d0, &
0.022838332222257d0, 0.695836086787803d0, 0.025734050548330d0, 0.116251915907597d0, &
0.116251915907597d0, 0.858014033544073d0, 0.858014033544073d0, 0.025734050548330d0, &
0.116251915907597d0, 0.025734050548330d0, 0.858014033544073d0, 0.116251915907597d0, &
0.025734050548330d0, 0.858014033544073d0 &
/)
w = (/ &
0.012865533220227d0, 0.012865533220227d0, 0.012865533220227d0, 0.021846272269019d0, &
0.021846272269019d0, 0.021846272269019d0, 0.031429112108943d0, 0.031429112108943d0, &
0.031429112108943d0, 0.017398056465355d0, 0.017398056465355d0, 0.017398056465355d0, &
0.003083130525780d0, 0.003083130525780d0, 0.003083130525780d0, 0.020185778883191d0, &
0.020185778883191d0, 0.020185778883191d0, 0.020185778883191d0, 0.020185778883191d0, &
0.020185778883191d0, 0.011178386601152d0, 0.011178386601152d0, 0.011178386601152d0, &
0.011178386601152d0, 0.011178386601152d0, 0.011178386601152d0, 0.008658115554329d0, &
0.008658115554329d0, 0.008658115554329d0, 0.008658115554329d0, 0.008658115554329d0, &
0.008658115554329d0 &
/)
elseif (p == 13) then
! Order 13d0 Dunavant Points
x = (/ &
0.333333333333333d0, 0.333333333333333d0, 0.009903630120591d0, 0.495048184939705d0, &
0.495048184939705d0, 0.495048184939705d0, 0.495048184939705d0, 0.009903630120591d0, &
0.062566729780852d0, 0.468716635109574d0, 0.468716635109574d0, 0.468716635109574d0, &
0.468716635109574d0, 0.062566729780852d0, 0.170957326397447d0, 0.414521336801277d0, &
0.414521336801277d0, 0.414521336801277d0, 0.414521336801277d0, 0.170957326397447d0, &
0.541200855914337d0, 0.229399572042831d0, 0.229399572042831d0, 0.229399572042831d0, &
0.229399572042831d0, 0.541200855914337d0, 0.771151009607340d0, 0.114424495196330d0, &
0.114424495196330d0, 0.114424495196330d0, 0.114424495196330d0, 0.771151009607340d0, &
0.950377217273082d0, 0.024811391363459d0, 0.024811391363459d0, 0.024811391363459d0, &
0.024811391363459d0, 0.950377217273082d0, 0.094853828379579d0, 0.268794997058761d0, &
0.268794997058761d0, 0.636351174561660d0, 0.636351174561660d0, 0.094853828379579d0, &
0.268794997058761d0, 0.094853828379579d0, 0.636351174561660d0, 0.268794997058761d0, &
0.094853828379579d0, 0.636351174561660d0, 0.018100773278807d0, 0.291730066734288d0, &
0.291730066734288d0, 0.690169159986905d0, 0.690169159986905d0, 0.018100773278807d0, &
0.291730066734288d0, 0.018100773278807d0, 0.690169159986905d0, 0.291730066734288d0, &
0.018100773278807d0, 0.690169159986905d0, 0.022233076674090d0, 0.126357385491669d0, &
0.126357385491669d0, 0.851409537834241d0, 0.851409537834241d0, 0.022233076674090d0, &
0.126357385491669d0, 0.022233076674090d0, 0.851409537834241d0, 0.126357385491669d0, &
0.022233076674090d0, 0.851409537834241d0 &
/)
w = (/ &
0.026260461700401d0, 0.005640072604665d0, 0.005640072604665d0, 0.005640072604665d0, &
0.015711759181227d0, 0.015711759181227d0, 0.015711759181227d0, 0.023536251252097d0, &
0.023536251252097d0, 0.023536251252097d0, 0.023681793268178d0, 0.023681793268178d0, &
0.023681793268178d0, 0.015583764522897d0, 0.015583764522897d0, 0.015583764522897d0, &
0.003987885732537d0, 0.003987885732537d0, 0.003987885732537d0, 0.018424201364366d0, &
0.018424201364366d0, 0.018424201364366d0, 0.018424201364366d0, 0.018424201364366d0, &
0.018424201364366d0, 0.008700731651911d0, 0.008700731651911d0, 0.008700731651911d0, &
0.008700731651911d0, 0.008700731651911d0, 0.008700731651911d0, 0.007760893419522d0, &
0.007760893419522d0, 0.007760893419522d0, 0.007760893419522d0, 0.007760893419522d0, &
0.007760893419522d0 &
/)
elseif (p == 14) then
! Order 14d0 Dunavant Points
x = (/ &
0.022072179275643d0, 0.488963910362179d0, 0.488963910362179d0, 0.488963910362179d0, &
0.488963910362179d0, 0.022072179275643d0, 0.164710561319092d0, 0.417644719340454d0, &
0.417644719340454d0, 0.417644719340454d0, 0.417644719340454d0, 0.164710561319092d0, &
0.453044943382323d0, 0.273477528308839d0, 0.273477528308839d0, 0.273477528308839d0, &
0.273477528308839d0, 0.453044943382323d0, 0.645588935174913d0, 0.177205532412543d0, &
0.177205532412543d0, 0.177205532412543d0, 0.177205532412543d0, 0.645588935174913d0, &
0.876400233818255d0, 0.061799883090873d0, 0.061799883090873d0, 0.061799883090873d0, &
0.061799883090873d0, 0.876400233818255d0, 0.961218077502598d0, 0.019390961248701d0, &
0.019390961248701d0, 0.019390961248701d0, 0.019390961248701d0, 0.961218077502598d0, &
0.057124757403648d0, 0.172266687821356d0, 0.172266687821356d0, 0.770608554774996d0, &
0.770608554774996d0, 0.057124757403648d0, 0.172266687821356d0, 0.057124757403648d0, &
0.770608554774996d0, 0.172266687821356d0, 0.057124757403648d0, 0.770608554774996d0, &
0.092916249356972d0, 0.336861459796345d0, 0.336861459796345d0, 0.570222290846683d0, &
0.570222290846683d0, 0.092916249356972d0, 0.336861459796345d0, 0.092916249356972d0, &
0.570222290846683d0, 0.336861459796345d0, 0.092916249356972d0, 0.570222290846683d0, &
0.014646950055654d0, 0.298372882136258d0, 0.298372882136258d0, 0.686980167808088d0, &
0.686980167808088d0, 0.014646950055654d0, 0.298372882136258d0, 0.014646950055654d0, &
0.686980167808088d0, 0.298372882136258d0, 0.014646950055654d0, 0.686980167808088d0, &
0.001268330932872d0, 0.118974497696957d0, 0.118974497696957d0, 0.879757171370171d0, &
0.879757171370171d0, 0.001268330932872d0, 0.118974497696957d0, 0.001268330932872d0, &
0.879757171370171d0, 0.118974497696957d0, 0.001268330932872d0, 0.879757171370171d0 &
/)
w = (/ &
0.010941790684715d0, 0.010941790684715d0, 0.010941790684715d0, 0.016394176772063d0, &
0.016394176772063d0, 0.016394176772063d0, 0.025887052253646d0, 0.025887052253646d0, &
0.025887052253646d0, 0.021081294368497d0, 0.021081294368497d0, 0.021081294368497d0, &
0.007216849834889d0, 0.007216849834889d0, 0.007216849834889d0, 0.002461701801200d0, &
0.002461701801200d0, 0.002461701801200d0, 0.012332876606282d0, 0.012332876606282d0, &
0.012332876606282d0, 0.012332876606282d0, 0.012332876606282d0, 0.012332876606282d0, &
0.019285755393531d0, 0.019285755393531d0, 0.019285755393531d0, 0.019285755393531d0, &
0.019285755393531d0, 0.019285755393531d0, 0.007218154056767d0, 0.007218154056767d0, &
0.007218154056767d0, 0.007218154056767d0, 0.007218154056767d0, 0.007218154056767d0, &
0.002505114419250d0, 0.002505114419250d0, 0.002505114419250d0, 0.002505114419250d0, &
0.002505114419250d0, 0.002505114419250d0 &
/)
elseif (p <= 17) then
! Order 17d0 Dunavant Points
x = (/ &
0.333333333333333d0, 0.333333333333333d0, 0.005658918886452d0, 0.497170540556774d0, &
0.497170540556774d0, 0.497170540556774d0, 0.497170540556774d0, 0.005658918886452d0, &
0.035647354750751d0, 0.482176322624625d0, 0.482176322624625d0, 0.482176322624625d0, &
0.482176322624625d0, 0.035647354750751d0, 0.099520061958437d0, 0.450239969020782d0, &
0.450239969020782d0, 0.450239969020782d0, 0.450239969020782d0, 0.099520061958437d0, &
0.199467521245206d0, 0.400266239377397d0, 0.400266239377397d0, 0.400266239377397d0, &
0.400266239377397d0, 0.199467521245206d0, 0.495717464058095d0, 0.252141267970953d0, &
0.252141267970953d0, 0.252141267970953d0, 0.252141267970953d0, 0.495717464058095d0, &
0.675905990683077d0, 0.162047004658461d0, 0.162047004658461d0, 0.162047004658461d0, &
0.162047004658461d0, 0.675905990683077d0, 0.848248235478508d0, 0.075875882260746d0, &
0.075875882260746d0, 0.075875882260746d0, 0.075875882260746d0, 0.848248235478508d0, &
0.968690546064356d0, 0.015654726967822d0, 0.015654726967822d0, 0.015654726967822d0, &
0.015654726967822d0, 0.968690546064356d0, 0.010186928826919d0, 0.334319867363658d0, &
0.334319867363658d0, 0.655493203809423d0, 0.655493203809423d0, 0.010186928826919d0, &
0.334319867363658d0, 0.010186928826919d0, 0.655493203809423d0, 0.334319867363658d0, &
0.010186928826919d0, 0.655493203809423d0, 0.135440871671036d0, 0.292221537796944d0, &
0.292221537796944d0, 0.572337590532020d0, 0.572337590532020d0, 0.135440871671036d0, &
0.292221537796944d0, 0.135440871671036d0, 0.572337590532020d0, 0.292221537796944d0, &
0.135440871671036d0, 0.572337590532020d0, 0.054423924290583d0, 0.319574885423190d0, &
0.319574885423190d0, 0.626001190286228d0, 0.626001190286228d0, 0.054423924290583d0, &
0.319574885423190d0, 0.054423924290583d0, 0.626001190286228d0, 0.319574885423190d0, &
0.054423924290583d0, 0.626001190286228d0, 0.012868560833637d0, 0.190704224192292d0, &
0.190704224192292d0, 0.796427214974071d0, 0.796427214974071d0, 0.012868560833637d0, &
0.190704224192292d0, 0.012868560833637d0, 0.796427214974071d0, 0.190704224192292d0, &
0.012868560833637d0, 0.796427214974071d0, 0.067165782413524d0, 0.180483211648746d0, &
0.180483211648746d0, 0.752351005937729d0, 0.752351005937729d0, 0.067165782413524d0, &
0.180483211648746d0, 0.067165782413524d0, 0.752351005937729d0, 0.180483211648746d0, &
0.067165782413524d0, 0.752351005937729d0, 0.014663182224828d0, 0.080711313679564d0, &
0.080711313679564d0, 0.904625504095608d0, 0.904625504095608d0, 0.014663182224828d0, &
0.080711313679564d0, 0.014663182224828d0, 0.904625504095608d0, 0.080711313679564d0, &
0.014663182224828d0, 0.904625504095608d0 &
/)
w = (/ &
0.016718599645402d0, 0.002546707720254d0, 0.002546707720254d0, 0.002546707720254d0, &
0.007335432263819d0, 0.007335432263819d0, 0.007335432263819d0, 0.012175439176836d0, &
0.012175439176836d0, 0.012175439176836d0, 0.015553775434484d0, 0.015553775434484d0, &
0.015553775434484d0, 0.015628555609310d0, 0.015628555609310d0, 0.015628555609310d0, &
0.012407827169833d0, 0.012407827169833d0, 0.012407827169833d0, 0.007028036535279d0, &
0.007028036535279d0, 0.007028036535279d0, 0.001597338086889d0, 0.001597338086889d0, &
0.001597338086889d0, 0.004059827659497d0, 0.004059827659497d0, 0.004059827659497d0, &
0.004059827659497d0, 0.004059827659497d0, 0.004059827659497d0, 0.013402871141582d0, &
0.013402871141582d0, 0.013402871141582d0, 0.013402871141582d0, 0.013402871141582d0, &
0.013402871141582d0, 0.009229996605411d0, 0.009229996605411d0, 0.009229996605411d0, &
0.009229996605411d0, 0.009229996605411d0, 0.009229996605411d0, 0.004238434267164d0, &
0.004238434267164d0, 0.004238434267164d0, 0.004238434267164d0, 0.004238434267164d0, &
0.004238434267164d0, 0.009146398385012d0, 0.009146398385012d0, 0.009146398385012d0, &
0.009146398385012d0, 0.009146398385012d0, 0.009146398385012d0, 0.003332816002083d0, &
0.003332816002083d0, 0.003332816002083d0, 0.003332816002083d0, 0.003332816002083d0, &
0.003332816002083d0 &
/)
elseif (p <= 19) then
! Order 19d0 Dunavant Points
x = (/ &
0.333333333333333d0, 0.333333333333333d0, 0.020780025853987d0, 0.489609987073006d0, &
0.489609987073006d0, 0.489609987073006d0, 0.489609987073006d0, 0.020780025853987d0, &
0.090926214604215d0, 0.454536892697893d0, 0.454536892697893d0, 0.454536892697893d0, &
0.454536892697893d0, 0.090926214604215d0, 0.197166638701138d0, 0.401416680649431d0, &
0.401416680649431d0, 0.401416680649431d0, 0.401416680649431d0, 0.197166638701138d0, &
0.488896691193805d0, 0.255551654403098d0, 0.255551654403098d0, 0.255551654403098d0, &
0.255551654403098d0, 0.488896691193805d0, 0.645844115695741d0, 0.177077942152130d0, &
0.177077942152130d0, 0.177077942152130d0, 0.177077942152130d0, 0.645844115695741d0, &
0.779877893544096d0, 0.110061053227952d0, 0.110061053227952d0, 0.110061053227952d0, &
0.110061053227952d0, 0.779877893544096d0, 0.888942751496321d0, 0.055528624251840d0, &
0.055528624251840d0, 0.055528624251840d0, 0.055528624251840d0, 0.888942751496321d0, &
0.974756272445543d0, 0.012621863777229d0, 0.012621863777229d0, 0.012621863777229d0, &
0.012621863777229d0, 0.974756272445543d0, 0.003611417848412d0, 0.395754787356943d0, &
0.395754787356943d0, 0.600633794794645d0, 0.600633794794645d0, 0.003611417848412d0, &
0.395754787356943d0, 0.003611417848412d0, 0.600633794794645d0, 0.395754787356943d0, &
0.003611417848412d0, 0.600633794794645d0, 0.134466754530780d0, 0.307929983880436d0, &
0.307929983880436d0, 0.557603261588784d0, 0.557603261588784d0, 0.134466754530780d0, &
0.307929983880436d0, 0.134466754530780d0, 0.557603261588784d0, 0.307929983880436d0, &
0.134466754530780d0, 0.557603261588784d0, 0.014446025776115d0, 0.264566948406520d0, &
0.264566948406520d0, 0.720987025817365d0, 0.720987025817365d0, 0.014446025776115d0, &
0.264566948406520d0, 0.014446025776115d0, 0.720987025817365d0, 0.264566948406520d0, &
0.014446025776115d0, 0.720987025817365d0, 0.046933578838178d0, 0.358539352205951d0, &
0.358539352205951d0, 0.594527068955871d0, 0.594527068955871d0, 0.046933578838178d0, &
0.358539352205951d0, 0.046933578838178d0, 0.594527068955871d0, 0.358539352205951d0, &
0.046933578838178d0, 0.594527068955871d0, 0.002861120350567d0, 0.157807405968595d0, &
0.157807405968595d0, 0.839331473680839d0, 0.839331473680839d0, 0.002861120350567d0, &
0.157807405968595d0, 0.002861120350567d0, 0.839331473680839d0, 0.157807405968595d0, &
0.002861120350567d0, 0.839331473680839d0, 0.223861424097916d0, 0.075050596975911d0, &
0.075050596975911d0, 0.701087978926173d0, 0.701087978926173d0, 0.223861424097916d0, &
0.075050596975911d0, 0.223861424097916d0, 0.701087978926173d0, 0.075050596975911d0, &
0.223861424097916d0, 0.701087978926173d0, 0.034647074816760d0, 0.142421601113383d0, &
0.142421601113383d0, 0.822931324069857d0, 0.822931324069857d0, 0.034647074816760d0, &
0.142421601113383d0, 0.034647074816760d0, 0.822931324069857d0, 0.142421601113383d0, &
0.034647074816760d0, 0.822931324069857d0, 0.010161119296278d0, 0.065494628082938d0, &
0.065494628082938d0, 0.924344252620784d0, 0.924344252620784d0, 0.010161119296278d0, &
0.065494628082938d0, 0.010161119296278d0, 0.924344252620784d0, 0.065494628082938d0, &
0.010161119296278d0, 0.924344252620784d0 &
/)
w = (/ &
0.016453165694459d0, 0.005165365945636d0, 0.005165365945636d0, 0.005165365945636d0, &
0.011193623631508d0, 0.011193623631508d0, 0.011193623631508d0, 0.015133062934734d0, &
0.015133062934734d0, 0.015133062934734d0, 0.015245483901099d0, 0.015245483901099d0, &
0.015245483901099d0, 0.012079606370821d0, 0.012079606370821d0, 0.012079606370821d0, &
0.008025401793400d0, 0.008025401793400d0, 0.008025401793400d0, 0.004042290130892d0, &
0.004042290130892d0, 0.004042290130892d0, 0.001039681013742d0, 0.001039681013742d0, &
0.001039681013742d0, 0.001942438452491d0, 0.001942438452491d0, 0.001942438452491d0, &
0.001942438452491d0, 0.001942438452491d0, 0.001942438452491d0, 0.012787080306011d0, &
0.012787080306011d0, 0.012787080306011d0, 0.012787080306011d0, 0.012787080306011d0, &
0.012787080306011d0, 0.004440451786669d0, 0.004440451786669d0, 0.004440451786669d0, &
0.004440451786669d0, 0.004440451786669d0, 0.004440451786669d0, 0.008062273380866d0, &
0.008062273380866d0, 0.008062273380866d0, 0.008062273380866d0, 0.008062273380866d0, &
0.008062273380866d0, 0.001245970908745d0, 0.001245970908745d0, 0.001245970908745d0, &
0.001245970908745d0, 0.001245970908745d0, 0.001245970908745d0, 0.009121420059476d0, &
0.009121420059476d0, 0.009121420059476d0, 0.009121420059476d0, 0.009121420059476d0, &
0.009121420059476d0, 0.005129281868099d0, 0.005129281868099d0, 0.005129281868099d0, &
0.005129281868099d0, 0.005129281868099d0, 0.005129281868099d0, 0.001899964427651d0, &
0.001899964427651d0, 0.001899964427651d0, 0.001899964427651d0, 0.001899964427651d0, &
0.001899964427651d0 &
/)
endif
xy = reshape(x,(/N,2/),order=(/2,1/))
end subroutine Gauss2D