Mesh Generation e Refinement

Mesh Generation e Refinement

Dopo aver studiato vari tipi di metodi utilizzati nella letteratura per i processi di mesh generation si è cercato di implementarne uno efficiente, affidabile e rapido. L’algoritmo costruito è stato pensato è testato per essere applicato a problemi di simulazioni agli elementi finiti per strutture di macchine elettriche innovative. Per questo motivo l’algoritmo prodotto doveva rispettare vincoli e confini spaziali che saranno dettagliati all’interno del capitolo. Oltre al codice essenziale per produrre una sequenza tetraedrica il lavoro ha previsto anche la costruzione di una interfaccia grafica adatta a risolvere problemi generici tridimensionali agli elementi finiti.

Tutte le tecniche precedentemente elencate di mesh refinement sono state implementate. Data una mesh non rifinita l’utente che utilizza il software può scegliere una qualunque delle tecniche di refinement.

Il Software è stato implementato su piattaforma Windows. La visualizzazione 3D presente all’interno del programma, è stata costruita attraverso il linguaggio di programmazione OpenGL [21, 22]. La possibilità di testare in ogni momento l’evoluzione della mesh tridimensionale ha permesso di migliorare passo per passo tutte le sequenze di calcolo automatico.

Nel campo degli algoritmi di mesh generation non sottoposti a vincoli di costruzioni (unstructured mesh), negli ultimi venti anni l’Advancing­front, ed i metodi di Delaunay hanno catalizzato l’attenzione e l’interesse dei maggiori centri di ricerca nel mondo [1, 2, 3].

Questi metodi implementano il processo di mesh generation in fasi differenti. Il primo passo compiuto da questi processi consiste nell’inserire tutti i nodi primari all’interno di una triangolazione provvisoria. Per includere ogni nodo nella mesh si potrebbe ad esempio costruire un esaedro diviso in cinque tetraedri. Il parallelepipedo deve essere grande tanto da contenere tutti i nodi primari nel suo volume interno. Durante il processo di costruzione un nuovo nodo da inserire si può trovare:

– All’interno di un tetraedro

– Su una faccia

– Su uno dei sei spigoli

Sono già stati elencati i metodi utilizzati per inserire ogni nuovo nodo in base ad ogni tipo di posizionamento spaziale.

Tuttavia la mesh finale corretta, deve contenere ogni nodo primario ma deve allo stesso tempo rispettare i vincoli imposti dalla struttura geometrica. La mesh provvisoria costruita inserendo I nodi in generale non conterrà al suo interno tutti i lati e le facce della struttura geometrica. In letteratura sono riportati un numero notevole di metodi che sono focalizzati alla modifica della mesh provvisoria necessaria al raggiungimento di tutte le condizioni richieste.

La maggior parte dei criteri di modifica della mesh provvisoria sono basati su delle modifiche locali di un numero ristretto di tetraedri: lo scopo di queste modifiche consiste nel sopprimere degli elementi che attraversano delle superfici di confine. Per sbrogliare dei domini molto complicati può anche essere necessario abbandonare il criterio di Delaunay e accettare anche triangolazioni differenti.

Hazlewood [6] ha dimostrato che è possibile produrre triangolazioni che rispettino i vincoli della struttura attraverso una ricostruzione locale della mesh. Il metodo proposto prevede la localizzazione della superficie di separazione tra la mesh da modificare e quella che invece non necessità correzione. Tutto lo spazio all’interno dal dominio da ricostruire deve essere sgombrato dai tetraedri che non permettono il rispetto dei vincoli strutturali, e successivamente la mesh può essere rigenerata utilizzando il criterio di Delaunay.

Altri metodi prevedono l’aggiunta di altri nodi secondari detti nodi di Steiner: se ad esempio uno spigolo di un tetraedro interseca una faccia triangolare di una struttura si potrebbe aggiungere un nuovo nodo coincidente con il punto di intersezione e dividere lo spigolo in due spigoli allineati. I nodi di Steiner in generale possono presentare però uno svantaggio: l’eccessiva suddivisione dei tetraedri in elementi più piccoli e con coefficienti di qualità spesso molto scarsi può provocare un collasso di alcuni tetraedri i cui quattro vertici possono diventare complanari.

Tutta la letteratura prodotta per cercare di raggiungere una mesh che rappresenti davvero una discretizzazione strutturale del problema FEM che si intende risolvere fanno però sorgere il dubbio che la problematica sia affetta da un errore di fondo: forse sarebbe più sensato costruire una mesh che in ogni momento cerchi di rispettare tutti i vincoli e le condizioni strutturali invece di cercare di modificare una mesh provvisoria per rispettare tutte le condizioni richieste.

L’algoritmo di mesh generation implementato è una variante del metodo di Delaunay – Voronoi costruito con dei controlli della struttura di confine effettuati durante l’evoluzione della mesh.

Utilizzando questo criterio, alla fine dell’inserzione di tutti i nodi non sarà più necessario modificare la reticolazione tetraedrica.

DESCRIZIONE DEL METODO DI TRIANGOLAZIONE AUTOMATICA

La struttura geometrica mostrata nella Figura 41 è uno schema di un polo elettrico di una macchina a flusso assiale. All’interno del dominio geometrico sono localizzate le parti fondamentali della macchina. In Figura 42 le superfici colorate in blu appartengono ai domini del ferro di rotore e statore della macchina. Il magnete permanente montato sul settore rotorico è posizionato tra i due elementi ferromagnetici.

Per applicare il metodo FE ogni settore della struttura deve avere delle proprietà fisiche ben definite (proprietà del materiale, condizioni al contorno, input prodotti dal magnete permanente). Ogni dominio strutturale è confinato da una mesh superficiale costituita da triangoli. Ovviamente ogni triangolo è definito da tre lati. I vertici del triangolo coincidono con tre nodi primari della geometria.

I nodi ed i lati della struttura sono disegnati nel wireframe della Figura 43. Alla fine della mesh generation ogni nodo della struttura deve coincidere almeno con un vertice di un tetraedro. Questa condizione è necessaria ma non sufficiente per ottenere una mesh definita e corretta.

Figura – Struttura di un polo elettrico di una macchina a flusso assiale. La struttura costruita e da reticolare è costituita da 17 settori.

Con il colore azzurro viene indicato lo spazio occupato dal ferro di statore e di rotore.

La struttura è composta da 132 nodi primari; 538 lati e 668 facce triangolari.

Figura – Wireframe della struttura. I nodi ed i lati mostrati devono essere poi inseriti tutti nella reticolazione della mesh.

Un’altra importante condizione è la seguente: ogni lato della struttura deve coincidere con almeno un lato di un tetraedro. Se le due condizioni sui nodi e sui lati sono corrette la mesh è corretta e il processo di generazione può considerarsi terminato.

Nella bibliografia il problema della inserzione ed eliminazione delle facce triangolari è ampiamente discusso. Comunque se i punti di Steiner non sono inseriti nella triangolazione questo caso può essere trascurato. Per dimostrare questa asserzione si consideri un dominio della struttura contornato da una mesh superficiale costituita da elementi triangolari. Se ogni lato della superficie è inserito nella mesh ed ogni lato di ogni tetraedro non interseca la superficie allora anche tutti gli elementi triangolari della superficie appartengono alla mesh 3D.

IL METODO PROPOSTO

L’algoritmo di mesh generation proposto prende spunto dal criterio di triangolazione Advancing-Front. Nella Figura 44 è riportato lo schema di un caso 2D. Il nodo della Figura 44 non è stato ancora inserito e si trova all’esterno della mesh in evoluzione il cui contorno esterno è visibile da una spezzata con spessore marcato. I quattro lati 1, 2, 3, 4 della mesh provvisoria sono affacciati al nodo da aggiungere. I nuovi triangoli possono essere formati congiungendo il nodo ai quattro lati. Il nodo viene quindi inserito all’interno della mesh insieme a quattro nuovi triangoli. Il metodo implementato al calcolatore è basato su questa idea di base. Questa esplicazione può essere facilmente estesa al caso 3D.

Durante la costruzione della mesh la superficie provvisoria esterna deve rimanere il più possibile convessa.

Figura – Inserzione di un nuovo nodo all’interno di una mesh bidimensionale

Per rispettare questa condizione di convessità e di compattezza della mesh in ogni momento del suo sviluppo la fase iniziale del processo prevede un ordinamento meticoloso della lista dei lati (e quindi dei nodi) della struttura da inserire. L’ordinamento dei nodi segue due criteri. Il primo tipo di ordinamento segue un criterio basato sulle distanza dai punti da inserire al baricentro della mesh 3D in costruzione. La lista dei lati può altrimenti essere ordinata dalla più piccola alla più elevata distanza tra i lati ed il punto medio del primo lato afferente al nodo iniziale. Un esempio di questo tipi di ordinamento è mostrato per un caso bidimensionale nella Figura 45a.

Il secondo tipo di ordinamento dei lati è implementato in modo da massimizzare in ogni momento il numero di lati interconnessi. Il risultato del secondo criterio di ordinamento è mostrato in Figura 45b. La visualizzazione di un ordinamento di una struttura 3D segue lo stesso criterio.

A questo punto è possibile iniziare il processo di generazione della triangolazione tetraedrica. Per cominciare si inserisce il primo tetraedro. Il primo tetraedro deve contenere almeno il primo lato strutturale della lista ordinata in precedenza.

Figura 45a,b: Esempio di ordinamento della mesh. a) Primo criterio. b) Secondo criterio.

 

Il primo tetraedro mostrato in Figura 46 presenta le seguenti proprietà: non contiene al suo interno nodi primari; non interseca nessun lato strutturale; non interseca nessuna faccia triangolare della struttura.

Figura – Wireframe della struttura e vista del primo tetraedro.

 

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

!*************************************************************************
Autori: 
Paolo Sordi Ordine degli ingegneri di Roma N° 22979 - Sezione A
Leonardo Santini Ordine degli ingegneri di Roma N° 22757 - Sezione A
 
!*************************************************************************
 
 
SUBROUTINE TRIANGOLAZIONE_A_CIPOLLONE(i_incastro)
USE IDENTIFIERS
USE CONFIG_MESH
USE CONFIG_PARAMETRI
IMPLICIT NONE
TYPE (WIN_MESSAGE) :: MESSAGE
INTEGER :: ITYPE
INTEGER je,kk
INTEGER num_lati_inseriti
LOGICAL inserito
INTEGER num_inserimenti
LOGICAL presente
LOGICAL incastrato
INTEGER i_incastro
INTEGER ntetra_old
INTEGER :: numero_lati_inseriti_x_tetra,num_lati_contenuti
DOUBLE PRECISION :: qualita_peggiore,qt
CALL WDIALOGSELECT (IDD_DIALOG001)
CALL WDIALOGSHOW   (-1,-1,0,MODELESS)
CALL CALCOLO_DEL_FATTORE_DI_SCALA_PER_LE_PROIEZIONI_ORTOGONALI ()
 
 
1000 CONTINUE
numero_lati_inseriti_x_tetra = 4
 
  DO JE = 1, NLAT
     CALL WMESSAGEPEEK (ITYPE, MESSAGE)
     IF (ITYPE.eq.CLOSEREQUEST) STOP
     CALL VERIFICA_SE_UN_LATO_PRIMARIO_SI_TROVA_NELLA_MESH (LATI (je,1),LATI (je,2), inserito)
     IF (.not.inserito) THEN
         CALL WDIALOGPUTINTEGER (idf_integer8,numero_lati_inseriti_x_tetra)
 
         CALL DISEGNA_LE_PROIEZIONI_ORTOGONALI_X_CARCIOFO (0)
         CALL DISEGNA_LE_PROIEZIONI_ORTOGONALI_DEI_TETRAEDRI ()
         CALL DISEGNA_LE_PROIEZIONI_DEL_LATO_DA_INSERIRE (LATI (je,1), LATI (je,2))
         CALL WDIALOGSELECT (IDD_DIALOG001)
         CALL WDIALOGPUTINTEGER (idf_integer2,NLAT)
         CALL WDIALOGPUTINTEGER (idf_integer7,je)
 
         CALL WDIALOGPUTINTEGER (idf_integer1,NNOD)
         CALL WDIALOGPUTINTEGER (idf_integer4,LATI (je,1))
         CALL WDIALOGPUTINTEGER (idf_integer5,LATI (je,2))
         CALL WDIALOGPUTINTEGER (idf_integer3 ,NTETRA)
!         CALL WDIALOGPUTINTEGER (idf_integer5 ,je)
         CALL WDIALOGPUTREAL    (IDF_REAL1, REAL (X(LATI (je,1))))
         CALL WDIALOGPUTREAL    (IDF_REAL2, REAL (Y(LATI (je,1))))
         CALL WDIALOGPUTREAL    (IDF_REAL3, REAL (Z(LATI (je,1))))
         CALL WDIALOGPUTREAL    (IDF_REAL4, REAL (X(LATI (je,2))))
         CALL WDIALOGPUTREAL    (IDF_REAL5, REAL (Y(LATI (je,2))))
         CALL WDIALOGPUTREAL    (IDF_REAL6, REAL (Z(LATI (je,2))))
         ntetra_old = NTETRA
         num_inserimenti = 0
         CALL VERIFICA_SE_IL_NODO_E_PRESENTE_NELLA_MESH (LATI (je,1), presente)
         IF (.NOT.presente) THEN
             CALL INSERIMENTO_DEL_NODO_NELLA_MESH (LATI (je,1),0)
             num_inserimenti = num_inserimenti + 1
         END IF
         CALL VERIFICA_SE_IL_NODO_E_PRESENTE_NELLA_MESH (LATI (je,2), presente)
         IF (.NOT.presente) THEN
             CALL INSERIMENTO_DEL_NODO_NELLA_MESH (LATI (je,2),0)
             num_inserimenti = num_inserimenti + 1
         END IF
         IF (i_incastro.eq.1) THEN
             CALL CONTROLLA_SE_QUALCHE_LATO_E_RIMASTO_INCASTRATO (incastrato)
             IF (incastrato) NTETRA = ntetra_old
             CALL CONTROLLA_SE_QUALCHE_FACCIA_E_RIMASTA_INCASTRATA (incastrato)
             IF (incastrato) NTETRA = ntetra_old
         END IF
! Verifica che il lato da inserire sia stato immesso nella mesh modificata
         inserito = .false.
         CALL VERIFICA_SE_UN_LATO_PRIMARIO_STA_NELLA_MESH (LATI (je,1),LATI (je,2),inserito)
         IF (.not.inserito)  NTETRA = ntetra_old
 
         IF (NTETRA.gt.ntetra_old) THEN
             CALL CONTROLLA_QUANTI_LATI_CONTIENE_IL_TETRA (NTETRA,num_lati_contenuti)
             IF (num_lati_contenuti.LT.numero_lati_inseriti_x_tetra) THEN
        !         NTETRA = ntetra_old
             END IF
         END IF
         IF (NTETRA.gt.ntetra_old) GOTO 1000
         CALL VERIFICA_SE_I_LATI_SONO_NELLA_MESH (num_lati_inseriti)
!         CALL WDIALOGPUTINTEGER (idf_integer26,num_lati_inseriti)
      END IF
  END DO
CALL WDIALOGHIDE()
END SUBROUTINE
 
 
!*************************************************************************
!*************************************************************************
 
SUBROUTINE  CONTROLLA_QUANTI_LATI_CONTIENE_IL_TETRA (num_tetra,num_lati_contenuti)
USE CONFIG_MESH
IMPLICIT NONE
INTEGER i1,i2,i3,i4,i5,i6,I
INTEGER :: num_lati_contenuti,num_tetra
I1=0
I2=0
I3=0
I4=0
i5=0
i6=0
DO I=1,NLAT
   IF (NT(num_tetra,1).EQ.LATI(i,1).AND.NT(num_tetra,2).EQ.LATI(i,2)) I1 = 1
   IF (NT(num_tetra,1).EQ.LATI(i,2).AND.NT(num_tetra,2).EQ.LATI(i,1)) I1 = 1
 
   IF (NT(num_tetra,1).EQ.LATI(i,1).AND.NT(num_tetra,3).EQ.LATI(i,2)) I2 = 1
   IF (NT(num_tetra,1).EQ.LATI(i,2).AND.NT(num_tetra,3).EQ.LATI(i,1)) I2 = 1
 
   IF (NT(num_tetra,1).EQ.LATI(i,1).AND.NT(num_tetra,4).EQ.LATI(i,2)) I3 = 1
   IF (NT(num_tetra,1).EQ.LATI(i,2).AND.NT(num_tetra,4).EQ.LATI(i,1)) I3 = 1
 
   IF (NT(num_tetra,2).EQ.LATI(i,1).AND.NT(num_tetra,3).EQ.LATI(i,2)) I4 = 1
   IF (NT(num_tetra,2).EQ.LATI(i,2).AND.NT(num_tetra,3).EQ.LATI(i,1)) I4 = 1
 
   IF (NT(num_tetra,2).EQ.LATI(i,1).AND.NT(num_tetra,4).EQ.LATI(i,2)) I5 = 1
   IF (NT(num_tetra,2).EQ.LATI(i,2).AND.NT(num_tetra,4).EQ.LATI(i,1)) I5 = 1
 
   IF (NT(num_tetra,3).EQ.LATI(i,1).AND.NT(num_tetra,4).EQ.LATI(i,2)) I6 = 1
   IF (NT(num_tetra,3).EQ.LATI(i,2).AND.NT(num_tetra,4).EQ.LATI(i,1)) I6 = 1
 
END DO
 
num_lati_contenuti = i1 + i2 + i3 + i4 + i5 + i6
 
END SUBROUTINE
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE INSERIMENTO_DEL_NODO_NELLA_MESH (nodo,i_nodo_in_mesh)
USE WINTERACTER
USE IDENTIFIERS
USE CONFIG_MESH
IMPLICIT NONE
INTEGER i_nodo_in_mesh
INTEGER nodo
INTEGER i
INTEGER ind
WRITE (1000,*) '****************** nodo da inserire ******************',nodo
CALL TROVA_SUPERFICIE_TRIANGOLARE_DELLA_MESH ()
CALL CHIAMATA_DEL_CICLO_DI_ATTESA_PER_SPLIT_CIPOLLONE (0,0)
DO i=1,NSURF
    CALL CONTROLLA_SE_TRIANGOLO_E_NODO_POSSONO_FORMARE_UN_TETRAEDRO (SURFACE(i,1),SURFACE(i,2),SURFACE(i,3),NODO,i_nodo_in_mesh)
END DO
END SUBROUTINE INSERIMENTO_DEL_NODO_NELLA_MESH
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE CONTROLLA_SE_TRIANGOLO_E_NODO_POSSONO_FORMARE_UN_TETRAEDRO (IN1,IN2,IN3,NODO,i_nodo_in_mesh)
USE CONFIG_MESH
IMPLICIT NONE
INTEGER IN1,IN2,IN3
INTEGER nodo
LOGICAL unisci
LOGICAL complanare
LOGICAL contiene
LOGICAL interseca
INTEGER i_nodo_in_mesh
IF (in1.eq.63.or.in2.eq.63.or.in3.eq.63) WRITE (1000,*) '63  ', IN1,IN2,IN3 ,'        nsurf  ', nsurf  ,'        Ntetra  ', Ntetra
IF (in1.eq.62.or.in2.eq.62.or.in3.eq.62) WRITE (1000,*) '62  ', IN1,IN2,IN3 ,'        nsurf  ', nsurf  ,'        Ntetra  ', Ntetra
 
CALL CONTROLLO_SPAZIALE_SE_TRIANGOLO_PUO_COLLEGARSI_UN_NODO (IN1,IN2,IN3,nodo,unisci)
IF (unisci) WRITE (1000,*) 'Questo triangolo potrebbe essere formato: ', IN1,IN2,IN3
IF (unisci) CALL CONTROLLA_NUM_TETRAEDRI_COMUNI_A_FACCIA (IN1,IN2,IN3,unisci)
IF (unisci) CALL CONTROLLA_NUM_TETRAEDRI_COMUNI_A_FACCIA (nodo,IN2,IN3,unisci)
IF (unisci) CALL CONTROLLA_NUM_TETRAEDRI_COMUNI_A_FACCIA (nodo,IN1,IN3,unisci)
IF (unisci) CALL CONTROLLA_NUM_TETRAEDRI_COMUNI_A_FACCIA (nodo,IN1,IN2,unisci)
IF (unisci) WRITE (1000,*) 'Il triangolo è di superficie ', IN1,IN2,IN3
!IF (.not.unisci) CALL WMESSAGEBOX (0,0,1,'Le altre facce del tetraedro sono doppie','Errore nella mesh')
IF (unisci) THEN
    NTETRA = NTETRA + 1
    NT(NTETRA,1) = IN1
    NT(NTETRA,2) = IN2
    NT(NTETRA,3) = IN3
    NT(NTETRA,4) = nodo
    CALL CONTROLLA_CHE_VERTICI_TETRAEDRO_NON_SIANO_COMPLANARI (NTETRA,complanare)
    IF (.not.complanare) THEN
        CALL ORDINANO (ntetra)
        IF (unisci) WRITE (1000,*) 'Il tetraedro non è complanare '
        CALL CONTROLLA_CHE_IL_TETRAEDRO_CONTENGA_UN_NODO (NTETRA,contiene)
        IF (contiene) THEN
            WRITE (1000,*) '!!!!NO!!!! Il tetraedro non è complanare '
            NTETRA = NTETRA - 1
            RETURN
        END IF
        CALL CONTROLLA_CHE_IL_TETRAEDRO_NON_INTERSECHI_DEI_LATI (NTETRA,interseca)
        IF (interseca) THEN
            WRITE (1000,*) '!!!!NO!!!! Il tetraedro interseca un lato '
            NTETRA = NTETRA - 1
            RETURN
        END IF
        CALL CONTROLLA_CHE_IL_TETRAEDRO_NON_INTERSECHI_UNA_FACCIA (NTETRA,interseca)
        IF (interseca) THEN
            WRITE (1000,*) '!!!!NO!!!! Il tetraedro interseca una faccia '
            NTETRA = NTETRA - 1
            RETURN
        END IF
    ELSE
        WRITE (1000,*) '!!!!NO!!!! Il tetraedro è complanare '
        NTETRA = NTETRA - 1
        RETURN
    END IF
 
    WRITE (1000,*) 'Unione con triangolo eseguita: faccia unita con il nodo',  NT(NTETRA,1),NT(NTETRA,2), NT(NTETRA,3)
END IF
END SUBROUTINE
 
 
 
 
!*************************************************************************
!*************************************************************************
 
 
 
SUBROUTINE CONTROLLO_SPAZIALE_SE_TRIANGOLO_PUO_COLLEGARSI_UN_NODO (in1,in2,in3,nodo,unisci)
USE WINTERACTER
USE IDENTIFIERS
USE CONFIG_MESH
IMPLICIT NONE
TYPE (WIN_MESSAGE) :: MESSAGE
INTEGER :: ITYPE
INTEGER :: nodo
LOGICAL :: unisci
INTEGER :: in1,in2,in3
INTEGER :: i
INTEGER :: ind1,ind2,ind3,ind4
DOUBLE PRECISION :: x_1,y_1,z_1,x_2,y_2,z_2,x_3,y_3,z_3,x_4,y_4,z_4
unisci = .true.
DO i=NTETRA,1,-1
   CALL WMESSAGEPEEK (ITYPE, MESSAGE)
   IF (ITYPE.eq.CLOSEREQUEST) STOP
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (NT (i,2),NT (i,4),NT (i,3),nodo,in1,x_1,y_1,z_1,IND1)
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (NT (i,1),NT (i,3),NT (i,4),nodo,in1,x_2,y_2,z_2,IND2)
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (NT (i,1),NT (i,4),NT (i,2),nodo,in1,x_3,y_3,z_3,IND3)
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (NT (i,1),NT (i,2),NT (i,3),nodo,in1,x_4,y_4,z_4,IND4)
   IF (ind1.ge.0.or.ind2.ge.0.or.ind3.ge.0.or.ind4.ge.0) THEN
      unisci = .false.
      RETURN
   END IF
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (NT (i,2),NT (i,4),NT (i,3),nodo,in2,x_1,y_1,z_1,IND1)
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (NT (i,1),NT (i,3),NT (i,4),nodo,in2,x_2,y_2,z_2,IND2)
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (NT (i,1),NT (i,4),NT (i,2),nodo,in2,x_3,y_3,z_3,IND3)
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (NT (i,1),NT (i,2),NT (i,3),nodo,in2,x_4,y_4,z_4,IND4)
   IF (ind1.ge.0.or.ind2.ge.0.or.ind3.ge.0.or.ind4.ge.0) THEN
      unisci = .false.
      RETURN
   END IF
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (NT (i,2),NT (i,4),NT (i,3),nodo,in3,x_1,y_1,z_1,IND1)
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (NT (i,1),NT (i,3),NT (i,4),nodo,in3,x_2,y_2,z_2,IND2)
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (NT (i,1),NT (i,4),NT (i,2),nodo,in3,x_3,y_3,z_3,IND3)
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (NT (i,1),NT (i,2),NT (i,3),nodo,in3,x_4,y_4,z_4,IND4)
   IF (ind1.ge.0.or.ind2.ge.0.or.ind3.ge.0.or.ind4.ge.0) THEN
      unisci = .false.
      RETURN
   END IF
END DO
END SUBROUTINE
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE CONTROLLA_CHE_IL_TETRAEDRO_NON_INTERSECHI_UNA_FACCIA (i_tetra,interseca)
USE CONFIG_MESH
IMPLICIT NONE
INTEGER i_tetra
LOGICAL interseca
INTEGER i
INTEGER ind1,ind2,ind3,ind4
INTEGER in1,in2,in3
DOUBLE PRECISION :: x_1,y_1,z_1
interseca = .false.
DO i=1,NFAC
   IN1 = FACCE (I,1)
   IN2 = FACCE (I,2)
   IN3 = FACCE (I,3)
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (IN1,IN2,IN3,NT (i_tetra,1),NT (i_tetra,2),x_1,y_1,z_1,IND1)
   IF (ind1.ge.0) THEN ; interseca = .true. ; RETURN ; END IF
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (IN1,IN2,IN3,NT (i_tetra,1),NT (i_tetra,3),x_1,y_1,z_1,IND1)
   IF (ind1.ge.0) THEN ; interseca = .true. ; RETURN ; END IF
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (IN1,IN2,IN3,NT (i_tetra,2),NT (i_tetra,3),x_1,y_1,z_1,IND1)
   IF (ind1.ge.0) THEN ; interseca = .true. ; RETURN ; END IF
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (IN1,IN2,IN3,NT (i_tetra,3),NT (i_tetra,4),x_1,y_1,z_1,IND1)
   IF (ind1.ge.0) THEN ; interseca = .true. ; RETURN ; END IF
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (IN1,IN2,IN3,NT (i_tetra,1),NT (i_tetra,4),x_1,y_1,z_1,IND1)
   IF (ind1.ge.0) THEN ; interseca = .true. ; RETURN ; END IF
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (IN1,IN2,IN3,NT (i_tetra,2),NT (i_tetra,4),x_1,y_1,z_1,IND1)
   IF (ind1.ge.0) THEN ; interseca = .true. ; RETURN ; END IF
END DO
END SUBROUTINE
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE  CONTROLLA_CHE_IL_TETRAEDRO_NON_INTERSECHI_DEI_LATI (i_tetra,interseca)
USE CONFIG_MESH
IMPLICIT NONE
INTEGER i_tetra
LOGICAL interseca
INTEGER i
INTEGER ind1,ind2,ind3,ind4
integer l1,l2
DOUBLE PRECISION x_1,y_1,z_1,x_2,y_2,z_2,x_3,y_3,z_3,x_4,y_4,z_4
interseca = .false.
DO i=1,NLAT
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (NT (i_tetra,2),NT (i_tetra,4),NT (i_tetra,3),LATI(i,1),LATI(i,2),x_1,y_1,z_1,IND1)
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (NT (i_tetra,1),NT (i_tetra,3),NT (i_tetra,4),LATI(i,1),LATI(i,2),x_2,y_2,z_2,IND2)
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (NT (i_tetra,1),NT (i_tetra,4),NT (i_tetra,2),LATI(i,1),LATI(i,2),x_3,y_3,z_3,IND3)
   CALL INCROCIO_TRA_TRIANGOLO_E_SEGMENTO (NT (i_tetra,1),NT (i_tetra,2),NT (i_tetra,3),LATI(i,1),LATI(i,2),x_4,y_4,z_4,IND4)
   l1 = LATI(i,1)
   l2 = LATI(i,2)
   IF (L1.EQ.NT(i_tetra,1).OR.L1.EQ.NT(i_tetra,2).OR.L1.EQ.NT(i_tetra,3).OR.L1.EQ.NT(i_tetra,4)) THEN
       IF (L2.EQ.NT(i_tetra,1).OR.L2.EQ.NT(i_tetra,2).OR.L2.EQ.NT(i_tetra,3).OR.L2.EQ.NT(i_tetra,4)) THEN
          IF (ind1.ge.0.or.ind2.ge.0.or.ind3.ge.0.or.ind4.ge.0) THEN
              CALL WMESSAGEBOX (0,0,1,'CONTROLLA_CHE_IL_TETRAEDRO_NON_INTERSECHI_DEI_LATI','Errore nella procedura')
          END IF
       END IF
   END IF
   IF (ind1.ge.0.or.ind2.ge.0.or.ind3.ge.0.or.ind4.ge.0) THEN
      interseca = .true.
      RETURN
   END IF
END DO
END SUBROUTINE
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE CONTROLLA_CHE_IL_TETRAEDRO_CONTENGA_UN_NODO (i_tetra,contiene)
USE CONFIG_MESH
IMPLICIT NONE
INTEGER I
LOGICAL contiene
INTEGER i_tetra
INTEGER IND
contiene = .FALSE.
DO i = 1 , NNOD
   CALL POSIZIONE_PUNTO_IE_RISPETTO_TETRAEDRO (i_tetra, x(I), y(I), z(I), ind, DBLE(1.E-2))
   IF (ind.GE.0) THEN
       contiene = .TRUE.
   END IF
END DO
END SUBROUTINE
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE  CONTROLLA_NUM_TETRAEDRI_COMUNI_A_FACCIA (in1,in2,in3,unisci)
USE CONFIG_MESH
IMPLICIT NONE
INTEGER i
INTEGER in1,in2,in3
INTEGER I_CONT
LOGICAL unisci
unisci = .true.
I_CONT = 0
DO i=1,NTETRA
   IF (in1.eq.NT(i,1).or.in1.eq.NT(i,2).or.in1.eq.NT(i,3).or.in1.eq.NT(i,4)) THEN
      IF (in2.eq.NT(i,1).or.in2.eq.NT(i,2).or.in2.eq.NT(i,3).or.in2.eq.NT(i,4)) THEN
          IF (in3.eq.NT(i,1).or.in3.eq.NT(i,2).or.in3.eq.NT(i,3).or.in3.eq.NT(i,4)) THEN
               I_CONT = I_CONT + 1
          END IF
      END IF
   END IF
END DO
IF (I_CONT.gt.1) unisci = .FALSE.
END SUBROUTINE
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE  TROVA_SUPERFICIE_TRIANGOLARE_DELLA_MESH ()
USE CONFIG_MESH
USE WINTERACTER
IMPLICIT NONE
LOGICAL unico
INTEGER I
nsurf = 0
!CALL WMESSAGEBOX (0,0,1,'INIZIO','')
DO i=1,NTETRA
     CALL CONTROLLA_NUM_TETRAEDRI_COMUNI_A_FACCIA (NT(i,1),NT(i,2),NT(i,3),unico)
     IF (unico) THEN
          nsurf = nsurf + 1
          surface (nsurf,1) = NT(i,2)
          surface (nsurf,2) = NT(i,1)
          surface (nsurf,3) = NT(i,3)
     END IF
     CALL CONTROLLA_NUM_TETRAEDRI_COMUNI_A_FACCIA (NT(i,1),NT(i,3),NT(i,4),unico)
     IF (unico) THEN
          nsurf = nsurf + 1
          surface (nsurf,1) = NT(i,3)
          surface (nsurf,2) = NT(i,1)
          surface (nsurf,3) = NT(i,4)
     END IF
     CALL CONTROLLA_NUM_TETRAEDRI_COMUNI_A_FACCIA (NT(i,1),NT(i,2),NT(i,4),unico)
     IF (unico) THEN
          nsurf = nsurf + 1
          surface (nsurf,1) = NT(i,4)
          surface (nsurf,2) = NT(i,1)
          surface (nsurf,3) = NT(i,2)
     END IF
     CALL CONTROLLA_NUM_TETRAEDRI_COMUNI_A_FACCIA (NT(i,2),NT(i,3),NT(i,4),unico)
     IF (unico) THEN
          nsurf = nsurf + 1
          surface (nsurf,1) = NT(i,4)
          surface (nsurf,2) = NT(i,2)
          surface (nsurf,3) = NT(i,3)
     END IF
END DO
!CALL WMESSAGEBOX (0,0,1,'FINE ','')
END SUBROUTINE
 
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE COSTRUISCI_IL_PRIMO_TETRAEDRO_DELLA_MESH ()
 
USE IDENTIFIERS
USE CONFIG_MESH
IMPLICIT NONE
LOGICAL :: complanare
LOGICAL :: contiene,interseca,interseca1
INTEGER :: i
CALL WDIALOGSELECT (IDD_DIALOG001)
CALL WDIALOGSHOW   (-1,-1,0,MODELESS)
CALL CALCOLO_DEL_FATTORE_DI_SCALA_PER_LE_PROIEZIONI_ORTOGONALI ()
 
NTETRA = 1
 
DO i=1,NFAC
   IF (FACCE(i,1).eq.1.or.FACCE(i,2).eq.1.or.FACCE(i,3).eq.1) THEN
       NT (1,1) = FACCE (i,1)
       NT (1,2) = FACCE (i,2)
       NT (1,3) = FACCE (i,3)
       GOTO 1000
   END IF
END  DO
 
!NT (1,1) = FACCE (1,1)
!NT (1,2) = FACCE (1,2)
!NT (1,3) = FACCE (1,3)
 
1000 CONTINUE
DO i=1,NLAT
    WRITE (1000,*) ' I vertici del 2° lato sono ',LATI (i,1),LATI (i,2)
    IF (LATI(i,2).eq.nt(1,1).or.LATI(i,2).eq.nt(1,2).or.LATI (i,2).eq.nt(1,3)) THEN
        NT (1,4) = LATI(i,1)
        CALL CONTROLLA_CHE_VERTICI_TETRAEDRO_NON_SIANO_COMPLANARI (1,complanare)
        IF (.not.complanare) THEN
            CALL ORDINANO (1)
            CALL CONTROLLA_CHE_IL_TETRAEDRO_CONTENGA_UN_NODO (1,contiene)
            CALL CONTROLLA_CHE_IL_TETRAEDRO_NON_INTERSECHI_DEI_LATI (1,interseca)
            CALL CONTROLLA_CHE_IL_TETRAEDRO_NON_INTERSECHI_UNA_FACCIA (1,interseca1)
            IF (.not.complanare.AND.(.NOT.CONTIENE).AND.(.NOT.INTERSECA).AND.(.NOT.INTERSECA1)) THEN
               WRITE (1000,*) ' Il quarto vertice del tetraedro è ',NT (1,4)
               CALL  DISEGNA_LE_PROIEZIONI_ORTOGONALI_X_CARCIOFO (0)
               CALL  DISEGNA_LE_PROIEZIONI_ORTOGONALI_DEI_TETRAEDRI ()
               CALL WDIALOGHIDE ()
               RETURN
            END IF
        END IF
    END IF
    IF (LATI(i,1).eq.nt(1,1).or.LATI(i,1).eq.nt(1,2).or.LATI (i,1).eq.nt(1,3)) THEN
        NT (1,4) = LATI(i,2)
        CALL CONTROLLA_CHE_VERTICI_TETRAEDRO_NON_SIANO_COMPLANARI (1,complanare)
        IF (.not.complanare) THEN
            CALL ORDINANO (1)
            CALL CONTROLLA_CHE_IL_TETRAEDRO_CONTENGA_UN_NODO (1,contiene)
            CALL CONTROLLA_CHE_IL_TETRAEDRO_NON_INTERSECHI_DEI_LATI (1,interseca)
            CALL CONTROLLA_CHE_IL_TETRAEDRO_NON_INTERSECHI_UNA_FACCIA (1,interseca1)
            IF (.not.complanare.AND.(.NOT.CONTIENE).AND.(.NOT.INTERSECA).AND.(.NOT.INTERSECA1)) THEN
               WRITE (1000,*) ' Il quarto vertice del tetraedro è ',NT (1,4)
               CALL  DISEGNA_LE_PROIEZIONI_ORTOGONALI_X_CARCIOFO (0)
               CALL  DISEGNA_LE_PROIEZIONI_ORTOGONALI_DEI_TETRAEDRI ()
               CALL WDIALOGHIDE ()
               RETURN
            END IF
        END IF
    END IF
END DO
END SUBROUTINE  
 
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE CONTROLLA_SE_QUALCHE_FACCIA_E_RIMASTA_INCASTRATA (incastrato)
 
USE IDENTIFIERS
USE CONFIG_MESH
IMPLICIT NONE
LOGICAL inserito
LOGICAL incastrato
LOGICAL nodo1_presente,nodo2_presente,nodo3_presente
INTEGER i,j
INTEGER l1,l2,l3
incastrato =.false.
DO i=1,NFAC
    nodo1_presente = .false.
    nodo2_presente = .false.
    nodo3_presente = .false.
    l1 = facce(i,1)
    l2 = facce(i,2)
    l3 = facce(i,3)
    CALL VERIFICA_SE_UN_FACCIA_SI_TROVA_NELLA_MESH (Facce(i,1),facce(i,2),facce(i,3),inserito)
    IF (.not.inserito) THEN
        DO j=1,NTETRA
            IF (L1.EQ.NT(j,1).OR.L1.EQ.NT(j,2).OR.L1.EQ.NT(j,3).OR.L1.EQ.NT(j,4)) THEN
               nodo1_presente = .true.
            END IF
            IF (L2.EQ.NT(j,1).OR.L2.EQ.NT(j,2).OR.L2.EQ.NT(j,3).OR.L2.EQ.NT(j,4)) THEN
               nodo2_presente = .true.
            END IF
            IF (L3.EQ.NT(j,1).OR.L3.EQ.NT(j,2).OR.L3.EQ.NT(j,3).OR.L3.EQ.NT(j,4)) THEN
               nodo3_presente = .true.
            END IF
        END DO
        IF (nodo1_presente.and.nodo2_presente.and.nodo3_presente) THEN
            incastrato = .TRUE.
!           CALL WMESSAGEBOX (0,0,1,'Faccia Incatrata','')
            RETURN
        END IF
    END IF
END DO
END SUBROUTINE
 
 
 
!*************************************************************************
!*************************************************************************
 
 
 
SUBROUTINE CONTROLLA_SE_QUALCHE_LATO_E_RIMASTO_INCASTRATO (incastrato)
USE WINTERACTER
USE IDENTIFIERS
USE CONFIG_MESH
IMPLICIT NONE
LOGICAL inserito
LOGICAL incastrato
LOGICAL nodo1_presente,nodo2_presente
INTEGER i,j
INTEGER l1,l2
incastrato =.false.
DO i=1,NLAT
    nodo1_presente = .false.
    nodo2_presente = .false.
    CALL VERIFICA_SE_UN_LATO_PRIMARIO_SI_TROVA_NELLA_MESH (LATI(i,1),LATI(i,2),inserito)
    IF (.not.inserito) THEN
        DO j=1,NTETRA
            l1 = LATI(i,1)
            l2 = LATI(i,2)
            IF (L1.EQ.NT(j,1).OR.L1.EQ.NT(j,2).OR.L1.EQ.NT(j,3).OR.L1.EQ.NT(j,4)) THEN
               nodo1_presente = .true.
            END IF
            IF (L2.EQ.NT(j,1).OR.L2.EQ.NT(j,2).OR.L2.EQ.NT(j,3).OR.L2.EQ.NT(j,4)) THEN
               nodo2_presente = .true.
            END IF
        END DO
        IF (nodo1_presente.and.nodo2_presente) THEN
            incastrato = .TRUE.
!           CALL WMESSAGEBOX (0,0,1,'Lato Incatrato','')
            RETURN
        END IF
    END IF
END DO
END SUBROUTINE
 
 
 
 
!*************************************************************************
!*************************************************************************
 
 
 
SUBROUTINE CHIAMATA_DEL_CICLO_DI_ATTESA_PER_SPLIT_CIPOLLONE (i_type,i_materiale)
USE WINTERACTER
USE IDENTIFIERS
USE CONFIG_MESH
USE CONFIG_PARAMETRI
INTEGER ITYPE
TYPE (WIN_MESSAGE) :: MESSAGE
CALL WDIALOGSELECT  (IDD_DIALOG_CIPOLLONE)
CALL WDIALOGGETCHECKBOX (IDF_CHECK3, ICHECKED)
IF (ICHECKED.eq.1) CALL DISEGNO3D()
CALL WDIALOGGETCHECKBOX (IDF_CHECK2, ICHECKED)
IF (ICHECKED.eq.1) THEN
   FERMA_SE_LATO_NON_INSERITO = .TRUE.
ELSE
   FERMA_SE_LATO_NON_INSERITO = .FALSE.
END IF
CALL WDIALOGGETCHECKBOX (IDF_CHECK1, ICHECKED)
IF (ICHECKED.eq.1) RETURN
DO
     CALL WMESSAGE(ITYPE, MESSAGE)
     IF (ITYPE.eq.PUSHBUTTON) THEN
          IF (MESSAGE%value1.eq.IDOK) THEN
              IF (i_type.eq.1) THEN
                  CALL WDIALOGGETINTEGER (IDF_INTEGER6,i_materiale)
              END IF
              RETURN
          END IF
     END IF
     IF (ITYPE.eq.CLOSEREQUEST) STOP
END DO
END SUBROUTINE

La mesh desiderata viene costruita tutta attorno a questo primo tetraedro. L’algoritmo di inserzione di ciascun lato strutturale può essere schematizzato con il diagramma di flusso della Figura 47. Il processo viene dettagliato di seguito.

  1. La lista dei lati viene scandagliata cercando il primo lato che non è stato ancora inserito nella triangolazione.

  2. Dato che la lista è stata ordinata, sicuramente uno dei due nodi è stato già inserito nella mesh provvisoria. Il nodo non ancora introdotto deve essere congiunto alla mesh con uno o più tetraedri di collegamento.

  3. Il nodo non inserito deve essere connesso alla superficie triangolare che contorna la mesh provvisoria: per creare un nuovo tetraedro il nodo viene collegato con una superficie triangolare affacciata al nodo stesso (terzo caso di Delaunay-Voronoy ).

  4. A questo punto il programma controlla se questo nuovo tetraedro può essere adatto alla triangolazione che si sta costruendo. Il nuovo tetraedro deve essere infatti eliminato s interseca una faccia o un lato della struttura. Il tetraedro non va ancora bene se contiene un nodo primario ancora da inserire.

  5. Ogni faccia triangolare del contorno della mesh viene controllata grazie ad un loop che ripete le istruzioni descritte ai punti 3 e 4.

Dopo aver scandagliato tutta la superficie esterna uno o più tetraedri sono stati inseriti nella mesh 3D. Questi elementi hanno sicuramente un vertice in comune che coincide con il nodo primario che si desiderava inserire.

La nuova configurazione non può comunque bloccare l’inserzione di alcuni lati o facce rimasti ancora da inserire. Può accadere ad esempio che i due estremi di un lato siano ad un certo punto inseriti all’interno della mesh quando invece il lato stesso non vi appartiene. In questo caso l’algoritmo è in grado di comprendere che la nuova configurazione non è corretta; I tetraedri devono essere eliminati e l’algoritmo deve cominciare ad introdurre nuovi lati della lista.

Figura 47   Diagramma a blocchi dell’algoritmo di mesh generation

APPLICAZIONI

La finestra di dialogo della figura Figura 48 riporta delle informazioni del progetto. In un campo numerico viene riportato il tempo necessario per produrre 508 tetraedri (13.8 secondi con un processore da 1300 MHz). I fattore di qualità permettono di apprezzare il grado di regolarità della triangolazione tetraedrica. In Figura 48 vengono riportati i tre coefficienti di qualità della mesh.

L’evoluzione del processo di mesh generation descritta finora può essere osservata nelle Figure 4-11,12,13. Il wireframe della mesh finale è mostrato in Figura 4-14a. La mesh prodotta costituisce una ripartizione corretta di ogni dominio. Di solito questa mesh di base non è adeguata per la soluzione di problemi FEM, come esplicato in precedenza è possibile applicare degli algoritmi per ottenere un buon livello di refinement.

Il processo di refinement è realizzato attraverso l’inserzione di nuovi nodi secondari. La disposizione di questi nuovi nodi è scelta per non abbattere drasticamente i valori del fattore di qualità della mesh.

Le Figure 4-15,16,17,18 riportano un processo di mesh generation di un polo di macchina a flusso assiale ad 8 poli. Una struttura non convessa costituita da un disco cilindrico è stata testata: i risultati sono riportati nelle Figure 4-19,20,21,22,23,24.

Le proiezioni ortogonali della mesh rifinita sono mostrate nella finestra di dialogo della Figure 12. Il software permette di selezionare qualsivoglia proiezione ortogonale per facilitare il processo di input dei dati.

Figura – Finestra di dialogo contenente le informazioni principali del problema.

.

Figura 49   Evoluzione della Mesh generation: 77 lati su 538 sono stati inseriti; 52 tetraedri sono stati creati.

Generazione mesh tetraedri fase intermedia
Evoluzione della Mesh: Superficie esterna – 226 lati su 538 sono stati inseriti; 209 tetraedri sono stati creati.

Figura 410   Evoluzione della Mesh: 226 lati su 538 sono stati inseriti; 209 tetraedri sono stati creati.

Figura 411   Evoluzione della Mesh generation: 441 lati su 538 sono stati inseriti; 443 tetraedri sono stati creati.

Figura 412a,b: Vista del wireframe della mesh. a) Mesh di base prima del refinement: 508 tetraedri sono stati creati,

b) Mesh dopo una fase del processo di refinement: 4478 elementi sono stati creati attraverso l’inserzione di 870 nodi secondari.

Figura – Wireframe della struttura e vista del primo tetraedro applicato ad una struttura di una macchina a flusso assiale ad 8 poli

Figura 414   Evoluzione della Mesh generation: 106 lati su 538 sono stati inseriti

Figura 415   Evoluzione della Mesh generation: 314 lati su 538 sono stati inseriti

.

Figura 416   Vista del wireframe della mesh. a) Mesh di base prima del refinement: 508 tetraedri sono stati create

 

.           

Figura 417   Disco cilindrico costituito da 140 nodi

(Filename: Esagono bucato 140 nodi regolare lati ordinati.net)

      

Figura 418   Evoluzione della Mesh generation: 21 lati su 420 sono stati inseriti

      

Figura 419   Evoluzione della Mesh generation: 97 lati su 420 sono stati inseriti

     

Figura 420   Evoluzione della Mesh generation: 203 lati su 420 sono stati inseriti

 

Figura 421   Evoluzione della Mesh generation: 302 lati su 420 sono stati inseriti

  

Figura 422   Evoluzione della Mesh generation: 410 lati su 420 sono stati inseriti

Figura – Proiezioni ortogonali di una struttura a scatola contenente un circuito magnetico ed un magnete permanente in rosso.

Il wireframe della mesh di base presenta i contorni in blu

Figura – Vista del wireframe della mesh. Il numero di tetraedri prodotti equivale a 59136.

Si noti la riduzione dei coefficienti di qualità globale rispetto alla triangolazione iniziale della Figura precedente

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
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
1001
1002
1003
1004
1005
1006
1007
1008
1009
1010
1011
1012
1013
1014
1015
1016
1017
1018
1019
1020
1021
1022
1023
1024
1025
1026
1027
1028
1029
1030
1031
1032
1033
1034
1035
1036
1037
1038
1039
1040
1041
1042
1043
1044
1045
1046
1047
1048
1049
1050

!*************************************************************************
Autori: 
Paolo Sordi Ordine degli ingegneri di Roma N° 22979 - Sezione A
Leonardo Santini Ordine degli ingegneri di Roma N° 22757 - Sezione A
 
!*************************************************************************
refinement.f90
 
SUBROUTINE RIFINISCI_LA_MESH_BIDIMENSIONALE_CON_CRITERIO_FRATTALE
USE CONFIG_MESH_2D
USE WINTERACTER
USE IDENTIFIERS
IMPLICIT NONE
TYPE (win_message) :: MESSAGE
INTEGER :: ITYPE 
LOGICAL :: continua
INTEGER :: ind,numlati
INTEGER :: It,i
INTEGER :: Iv1,iv2,iv3,im1,im2,im3,im4,im5,im6
INTEGER, DIMENSION (NTRI_2D*3,3) :: elenco_lati
CALL WDIALOGPUTINTEGER (idf_integer2 , NTRI_2D*8)
 
CALL CALCOLO_DEL_FATTORE_DI_SCALA_PER_LE_PROIEZIONI_ORTOGONALI ()
numlati = 0
DO It=1, NTRI_2D
    CALL INSERISCI_UN_NODO_IN_MEZZO_AL_LATO_X_BIDIMENSIONALE (NT_2D(it,1),NT_2D(it,2), elenco_lati, numlati)
    CALL INSERISCI_UN_NODO_IN_MEZZO_AL_LATO_X_BIDIMENSIONALE (NT_2D(it,1),NT_2D(it,3), elenco_lati, numlati)
    CALL INSERISCI_UN_NODO_IN_MEZZO_AL_LATO_X_BIDIMENSIONALE (NT_2D(it,2),NT_2D(it,3), elenco_lati, numlati)
END DO
 
DO It=1, NTRI_2D
    iv1 =  NT_2D(it,1)
    iv2 =  NT_2D(it,2)
    iv3 =  NT_2D(it,3)
    DO i=1,numlati
        IF (elenco_lati (i,1).eq.iv1.and.elenco_lati (i,2).eq.iv2) im1 = elenco_lati (i,3)
        IF (elenco_lati (i,2).eq.iv1.and.elenco_lati (i,1).eq.iV2) im1 = elenco_lati (i,3)
 
        IF (elenco_lati (i,1).eq.iv1.and.elenco_lati (i,2).eq.iv3) im2 = elenco_lati (i,3)
        IF (elenco_lati (i,2).eq.iv1.and.elenco_lati (i,1).eq.iV3) im2 = elenco_lati (i,3)
 
        IF (elenco_lati (i,1).eq.iv2.and.elenco_lati (i,2).eq.iv3) im3 = elenco_lati (i,3)
        IF (elenco_lati (i,2).eq.iv2.and.elenco_lati (i,1).eq.iV3) im3 = elenco_lati (i,3)
        CALL WMESSAGEPEEK (ITYPE, MESSAGE)
        IF (ITYPE.eq.CLOSEREQUEST) STOP
    END DO
    CALL WMESSAGEPEEK (ITYPE, MESSAGE)
    IF (ITYPE.eq.CLOSEREQUEST) STOP
    NT_2D(it,1) =  iv1
    NT_2D(it,2) =  im1
    NT_2D(it,3) =  im2
    NTRI_2D = NTRI_2D + 1
    NT_2D(NTRI_2D,1) =  iv2
    NT_2D(NTRI_2D,2) =  im3
    NT_2D(NTRI_2D,3) =  im1
    NTRI_2D = NTRI_2D + 1
    NT_2D(NTRI_2D,1) =  iv3
    NT_2D(NTRI_2D,2) =  im2
    NT_2D(NTRI_2D,3) =  im3
    NTRI_2D = NTRI_2D + 1
    NT_2D(NTRI_2D,1) =  im1
    NT_2D(NTRI_2D,2) =  im3
    NT_2D(NTRI_2D,3) =  im2
!    CALL WDIALOGPUTINTEGER (idf_integer3 , NTETRA)
!    CALL WDIALOGPUTINTEGER (idf_integer2 , NLAT)
!    CALL WDIALOGPUTINTEGER (idf_integer1 , NNOD)
    CALL WMESSAGEPEEK      (ITYPE, MESSAGE)
    IF (ITYPE.eq.CLOSEREQUEST) STOP
END DO
END SUBROUTINE
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE INSERISCI_UN_NODO_IN_MEZZO_AL_LATO_X_BIDIMENSIONALE (in1, in2, elenco_lati, numlati)
USE CONFIG_MESH_2D
IMPLICIT NONE
INTEGER, DIMENSION (NTRI_2D*3,3) :: elenco_lati
INTEGER :: in1,in2,numlati
INTEGER :: i
DO i=1,numlati
    if (elenco_lati (i,1).eq.in1.and.elenco_lati (i,2).eq.in2) RETURN
    if (elenco_lati (i,2).eq.in1.and.elenco_lati (i,1).eq.in2) RETURN
END DO
NNOD_2D = NNOD_2D + 1
X_2D (NNOD_2D) = (X_2D(in1) + X_2D(in2)) / 2.0
Y_2D (NNOD_2D) = (Y_2D(in1) + Y_2D(in2)) / 2.0
Z_2D (NNOD_2D) = (z_2D(in1) + z_2D(in2)) / 2.0
numlati = numlati + 1
elenco_lati (numlati,1) = in1
elenco_lati (numlati,2) = in2
elenco_lati (numlati,3) = NNOD_2D
END SUBROUTINE
 
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE RIFINISCI_LA_MESH_CON_CRITERIO_FRATTALE (i_esci)
USE CONFIG_MESH
USE CONFIG_FEM
USE WINTERACTER
USE IDENTIFIERS
IMPLICIT NONE
TYPE (win_message) :: MESSAGE
INTEGER :: ITYPE,i_esci
LOGICAL :: continua
INTEGER :: ind,numlati
INTEGER :: It,i
INTEGER :: Iv1,iv2,iv3,iv4,im1,im2,im3,im4,im5,im6
 
 
INTEGER, DIMENSION (MAXTET*3,3) :: elenco_lati
IF (i_esci.ne.1) CALL WDIALOGSELECT     (IDD_DIALOG001)
IF (i_esci.ne.1) CALL WDIALOGPUTINTEGER (idf_integer5 , NTETRA*8)
 
CALL CALCOLO_DEL_FATTORE_DI_SCALA_PER_LE_PROIEZIONI_ORTOGONALI ()
numlati = 0
 
CALL WDIALOGRANGEPROGRESSBAR (IDF_PROGRESS1,1,NTETRA)
 
 
DO It = 1, NTETRA
    CALL WDIALOGPUTPROGRESSBAR (IDF_PROGRESS1,it,Absolute)
 
    CALL INSERISCI_UN_NODO_IN_MEZZO_AL_LATO (NT(it,1),NT(it,2), elenco_lati, numlati)
    CALL INSERISCI_UN_NODO_IN_MEZZO_AL_LATO (NT(it,1),NT(it,3), elenco_lati, numlati)
    CALL INSERISCI_UN_NODO_IN_MEZZO_AL_LATO (NT(it,1),NT(it,4), elenco_lati, numlati)
    CALL INSERISCI_UN_NODO_IN_MEZZO_AL_LATO (NT(it,2),NT(it,3), elenco_lati, numlati)
    CALL INSERISCI_UN_NODO_IN_MEZZO_AL_LATO (NT(it,2),NT(it,4), elenco_lati, numlati)
    CALL INSERISCI_UN_NODO_IN_MEZZO_AL_LATO (NT(it,3),NT(it,4), elenco_lati, numlati)
    CALL WDIALOGPUTINTEGER (idf_integer1, NNOD)
    CALL WMESSAGEPEEK      (ITYPE, MESSAGE)
    IF (ITYPE.eq.CLOSEREQUEST) STOP
END DO
 
IF (i_esci.eq.1) RETURN
 
DO It = 1, NTETRA
    iv1 =  NT(it,1)
    iv2 =  NT(it,2)
    iv3 =  NT(it,3)
    iv4 =  NT(it,4)
    DO i=1,numlati
        IF (elenco_lati (i,1).eq.iv1.and.elenco_lati (i,2).eq.iv2) im1 = elenco_lati (i,3)
        IF (elenco_lati (i,2).eq.iv1.and.elenco_lati (i,1).eq.iV2) im1 = elenco_lati (i,3)
 
        IF (elenco_lati (i,1).eq.iv1.and.elenco_lati (i,2).eq.iv3) im2 = elenco_lati (i,3)
        IF (elenco_lati (i,2).eq.iv1.and.elenco_lati (i,1).eq.iV3) im2 = elenco_lati (i,3)
 
        IF (elenco_lati (i,1).eq.iv1.and.elenco_lati (i,2).eq.iv4) im3 = elenco_lati (i,3)
        IF (elenco_lati (i,2).eq.iv1.and.elenco_lati (i,1).eq.iV4) im3 = elenco_lati (i,3)
 
        IF (elenco_lati (i,1).eq.iv2.and.elenco_lati (i,2).eq.iv3) im4 = elenco_lati (i,3)
        IF (elenco_lati (i,2).eq.iv2.and.elenco_lati (i,1).eq.iV3) im4 = elenco_lati (i,3)
 
        IF (elenco_lati (i,1).eq.iv2.and.elenco_lati (i,2).eq.iv4) im5 = elenco_lati (i,3)
        IF (elenco_lati (i,2).eq.iv2.and.elenco_lati (i,1).eq.iV4) im5 = elenco_lati (i,3)
 
        IF (elenco_lati (i,1).eq.iv3.and.elenco_lati (i,2).eq.iv4) im6 = elenco_lati (i,3)
        IF (elenco_lati (i,2).eq.iv3.and.elenco_lati (i,1).eq.iV4) im6 = elenco_lati (i,3)
        CALL WMESSAGEPEEK (ITYPE, MESSAGE)
        IF (ITYPE.eq.CLOSEREQUEST) STOP
    END DO
    CALL WMESSAGEPEEK (ITYPE, MESSAGE)
    IF (ITYPE.eq.CLOSEREQUEST) STOP
    NT(it,1) =  iv1
    NT(it,2) =  im1
    NT(it,3) =  im2
    NT(it,4) =  im3
 
    CALL ORDINANO (it)
    NTETRA = NTETRA + 1
    NT(NTETRA,1) =  iv2
    NT(NTETRA,2) =  im1
    NT(NTETRA,3) =  im4
    NT(NTETRA,4) =  im5
    NT_PRIMARIO (NTETRA) =  NT_PRIMARIO (it)
 
 
    CALL ORDINANO (NTETRA)
    MATERIALE_TETRAEDRO (NTETRA) = MATERIALE_TETRAEDRO (it)
    NTETRA = NTETRA + 1
    NT(NTETRA,1) =  iv3
    NT(NTETRA,2) =  im2
    NT(NTETRA,3) =  im4
    NT(NTETRA,4) =  im6
    NT_PRIMARIO (NTETRA) =  NT_PRIMARIO (it)
 
    CALL ORDINANO (NTETRA)
    MATERIALE_TETRAEDRO (NTETRA) = MATERIALE_TETRAEDRO (it)
    NTETRA = NTETRA + 1
    NT(NTETRA,1) =  iv4
    NT(NTETRA,2) =  im3
    NT(NTETRA,3) =  im5
    NT(NTETRA,4) =  im6
    NT_PRIMARIO (NTETRA) = NT_PRIMARIO (it)
    CALL ORDINANO (NTETRA)
    MATERIALE_TETRAEDRO (NTETRA) = MATERIALE_TETRAEDRO (it)
!-------------------------------------------------------------------------
    NTETRA = NTETRA + 1
    NT(NTETRA,1) =  im1
    NT(NTETRA,2) =  im6
    NT(NTETRA,3) =  im2
    NT(NTETRA,4) =  im3
    NT_PRIMARIO (NTETRA) = NT_PRIMARIO (it)
    CALL ORDINANO (NTETRA)
    MATERIALE_TETRAEDRO (NTETRA) = MATERIALE_TETRAEDRO (it)
    NTETRA = NTETRA + 1
    NT(NTETRA,1) =  im1
    NT(NTETRA,2) =  im6
    NT(NTETRA,3) =  im2
    NT(NTETRA,4) =  im4
    NT_PRIMARIO (NTETRA) = NT_PRIMARIO (it)
    CALL ORDINANO (NTETRA)
    MATERIALE_TETRAEDRO (NTETRA) = MATERIALE_TETRAEDRO (it)
    NTETRA = NTETRA + 1
    NT(NTETRA,1) =  im1
    NT(NTETRA,2) =  im6
    NT(NTETRA,3) =  im5
    NT(NTETRA,4) =  im4
    NT_PRIMARIO (NTETRA) = NT_PRIMARIO (it)
    CALL ORDINANO (NTETRA)
    MATERIALE_TETRAEDRO (NTETRA) = MATERIALE_TETRAEDRO (it)
    NTETRA = NTETRA + 1
    NT(NTETRA,1) =  im1
    NT(NTETRA,2) =  im6
    NT(NTETRA,3) =  im3
    NT(NTETRA,4) =  im5
    NT_PRIMARIO (NTETRA) = NT_PRIMARIO (it)
    CALL ORDINANO (NTETRA)
    MATERIALE_TETRAEDRO (NTETRA) = MATERIALE_TETRAEDRO (it)
    CALL WDIALOGPUTINTEGER (idf_integer3 , NTETRA)
    CALL WDIALOGPUTINTEGER (idf_integer2 , NLAT)
    CALL WDIALOGPUTINTEGER (idf_integer1 , NNOD)
    CALL WMESSAGEPEEK      (ITYPE, MESSAGE)
    IF (ITYPE.eq.CLOSEREQUEST) STOP
END DO
CALL DISEGNA_LE_PROIEZIONI_ORTOGONALI_X_CARCIOFO (0)
CALL DISEGNA_LE_PROIEZIONI_ORTOGONALI_DEI_TETRAEDRI ()
CALL WDIALOGHIDE
 
!call SCRITTURA_DELLE_COORDINATE_SU_FILE()
!pause
END SUBROUTINE
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE INSERISCI_UN_NODO_IN_MEZZO_AL_LATO (in1, in2, elenco_lati, numlati)
USE CONFIG_FEM
USE CONFIG_MESH
IMPLICIT NONE
INTEGER, DIMENSION (MAXTET*3,3) :: elenco_lati
INTEGER :: in1,in2,numlati
INTEGER :: i
DO i=1,numlati
    if (elenco_lati (i,1).eq.in1.and.elenco_lati (i,2).eq.in2) RETURN
    if (elenco_lati (i,2).eq.in1.and.elenco_lati (i,1).eq.in2) RETURN
END DO
NNOD = NNOD + 1
CONDIZIONE_INIZIALE (NNOD) = -1.0
X (NNOD) = (X(in1) + X(in2)) / 2.0
Y (NNOD) = (Y(in1) + Y(in2)) / 2.0
Z (NNOD) = (z(in1) + z(in2)) / 2.0
numlati = numlati + 1
elenco_lati (numlati,1) = in1
elenco_lati (numlati,2) = in2
elenco_lati (numlati,3) = nnod
END SUBROUTINE
 
 
 
!*************************************************************************
!*************************************************************************
 
 
 
SUBROUTINE SPLITTA_IL_TETRAEDRO_SUL_LATO_PIU_LUNGO ()
USE CONFIG_MESH
IMPLICIT NONE
INTEGER :: k_div
 
INTEGER :: k_migliore
INTEGER :: i_best,j_best,it_best
 
INTEGER :: i_tetra
DOUBLE PRECISION :: lunghezza_best
DOUBLE PRECISION :: x1,y1,z1,x2,y2,z2
 
DO i_tetra = 1,NTETRA
   lunghezza_best = 0.0
   i_best = 0
   j_best = 0
   CALL TROVA_LATO_PIU_LUNGO (i_tetra,1,2,lunghezza_best,i_best,j_best,it_best)
   CALL TROVA_LATO_PIU_LUNGO (i_tetra,1,3,lunghezza_best,i_best,j_best,it_best)
   CALL TROVA_LATO_PIU_LUNGO (i_tetra,1,4,lunghezza_best,i_best,j_best,it_best)
   CALL TROVA_LATO_PIU_LUNGO (i_tetra,2,3,lunghezza_best,i_best,j_best,it_best)
   CALL TROVA_LATO_PIU_LUNGO (i_tetra,2,4,lunghezza_best,i_best,j_best,it_best)
   CALL TROVA_LATO_PIU_LUNGO (i_tetra,3,4,lunghezza_best,i_best,j_best,it_best)
   x1 = X(NT(i_tetra,i_best))
   y1 = Y(NT(i_tetra,i_best))
   z1 = Z(NT(i_tetra,i_best))
   x2 = X(NT(i_tetra,j_best))
   y2 = Y(NT(i_tetra,j_best))
   z2 = Z(NT(i_tetra,j_best))
   NNOD = NNOD + 1
   X (NNOD) = x1 + (x2-x1) / 2.0
   Y (NNOD) = y1 + (y2-y1) / 2.0
   Z (NNOD) = z1 + (z2-z1) / 2.0
   CALL TETRAEDIZZAZIONE_PER_NODI (NNOD, 0, 0, NNOD)
END DO
CALL WDIALOGHIDE
END SUBROUTINE
 
 
 
 
 
!*************************************************************************
!*************************************************************************
 
 
 
SUBROUTINE TROVA_I_TETRAEDRI_CHE_CONTENGONO_LO_SPIGOLO (it_array,i_vert_1,i_vert_2,i_vert_opp1,i_vert_opp2,trovato)
USE CONFIG_MESH
USE CONFIG_PARAMETRI
USE IDENTIFIERS
USE WINTERACTER
IMPLICIT NONE
INTEGER,DIMENSION (0:100) :: it_array
LOGICAL :: trovato
INTEGER :: i
INTEGER :: i_vert_1,i_vert_2
INTEGER :: i_vert_opp1,i_vert_opp2
INTEGER :: num_tetraedri
 
TETRA_IN_STRUCTURE = .true.
 
trovato =.false.
it_array = 0
DO i=1,NTETRA
   IF (i_vert_1.EQ.NT(i,1).OR.i_vert_1.EQ.NT(i,2).OR.i_vert_1.EQ.NT(i,3).OR.i_vert_1.EQ.NT(i,4)) THEN
       IF (i_vert_2.EQ.NT(i,1).OR.i_vert_2.EQ.NT(i,2).OR.i_vert_2.EQ.NT(i,3).OR.i_vert_2.EQ.NT(i,4)) THEN
            TETRA_IN_STRUCTURE (i) = .false.
            trovato = .true.
            it_array(1) = i
            IF (NT(i,1).NE.i_vert_1.AND.NT(i,1).NE.i_vert_2) i_vert_opp1 = NT (i,1)
            IF (NT(i,2).NE.i_vert_1.AND.NT(i,2).NE.i_vert_2) i_vert_opp1 = NT (i,2)
            IF (NT(i,3).NE.i_vert_1.AND.NT(i,3).NE.i_vert_2) i_vert_opp1 = NT (i,3)
            IF (NT(i,4).NE.i_vert_1.AND.NT(i,4).NE.i_vert_2) i_vert_opp1 = NT (i,4)
            IF (NT(i,1).NE.i_vert_1.AND.NT(i,1).NE.i_vert_2.AND.NT(i,1).NE.i_vert_opp1) i_vert_opp2 = NT (i,1)
            IF (NT(i,2).NE.i_vert_1.AND.NT(i,2).NE.i_vert_2.AND.NT(i,2).NE.i_vert_opp1) i_vert_opp2 = NT (i,2)
            IF (NT(i,3).NE.i_vert_1.AND.NT(i,3).NE.i_vert_2.AND.NT(i,3).NE.i_vert_opp1) i_vert_opp2 = NT (i,3)
            IF (NT(i,4).NE.i_vert_1.AND.NT(i,4).NE.i_vert_2.AND.NT(i,4).NE.i_vert_opp1) i_vert_opp2 = NT (i,4)
            GOTO 100
       END IF
   END IF
END DO
100 CONTINUE
num_tetraedri = 1
DO i=1,NTETRA
   IF (i_vert_1.EQ.NT(i,1).OR.i_vert_1.EQ.NT(i,2).OR.i_vert_1.EQ.NT(i,3).OR.i_vert_1.EQ.NT(i,4)) THEN
       IF (i_vert_2.EQ.NT(i,1).OR.i_vert_2.EQ.NT(i,2).OR.i_vert_2.EQ.NT(i,3).OR.i_vert_2.EQ.NT(i,4)) THEN
           IF (i.ne.it_array(1)) THEN
                TETRA_IN_STRUCTURE (i) = .false.
                num_tetraedri = num_tetraedri + 1
                it_array(0) = num_tetraedri
                it_array(num_tetraedri) = i
           END IF
       END IF
   END IF
END DO
!IF (num_tetraedri.gt.3) THEN
!    write (1000,*) 'Numero di tetredri per lo spigolo .ne. 3',num_tetraedri
!    CALL WDIALOGPUTINTEGER (idf_integer37,num_tetraedri)
!END IF
END SUBROUTINE
 
 
 
!*************************************************************************
!*************************************************************************
 
 
 
SUBROUTINE TROVA_I_DUE_TETRAEDRI_CHE_CONTENGANO_LA_FACCIA (i_tetra1, i_tetra2, i_vert_1, i_vert_2, i_vert_3, i_vert_opp, trovato)
USE CONFIG_MESH
USE CONFIG_PARAMETRI
IMPLICIT NONE
LOGICAL :: trovato
INTEGER :: i
INTEGER :: i_tetra1,i_tetra2
INTEGER :: i_vert_1,i_vert_2,i_vert_3
INTEGER :: i_vert_opp
trovato =.false.
i_tetra1= 0
i_tetra2= 0
TETRA_IN_STRUCTURE = .true.
DO i=1,NTETRA
   IF (i_vert_1.EQ.NT(i,1).OR.i_vert_1.EQ.NT(i,2).OR.i_vert_1.EQ.NT(i,3).OR.i_vert_1.EQ.NT(i,4)) THEN
       IF (i_vert_2.EQ.NT(i,1).OR.i_vert_2.EQ.NT(i,2).OR.i_vert_2.EQ.NT(i,3).OR.i_vert_2.EQ.NT(i,4)) THEN
           IF (i_vert_3.EQ.NT(i,1).OR.i_vert_3.EQ.NT(i,2).OR.i_vert_3.EQ.NT(i,3).OR.i_vert_3.EQ.NT(i,4)) THEN
               trovato = .true.
               i_tetra1 = i
               TETRA_IN_STRUCTURE (i) = .false.
               IF (NT(i_tetra1,1).NE.i_vert_1.AND.NT(i_tetra1,1).NE.i_vert_2.AND.NT(i_tetra1,1).NE.i_vert_3)  i_vert_opp = NT (i_tetra1,1)
               IF (NT(i_tetra1,2).NE.i_vert_1.AND.NT(i_tetra1,2).NE.i_vert_2.AND.NT(i_tetra1,2).NE.i_vert_3)  i_vert_opp = NT (i_tetra1,2)
               IF (NT(i_tetra1,3).NE.i_vert_1.AND.NT(i_tetra1,3).NE.i_vert_2.AND.NT(i_tetra1,3).NE.i_vert_3)  i_vert_opp = NT (i_tetra1,3)
               IF (NT(i_tetra1,4).NE.i_vert_1.AND.NT(i_tetra1,4).NE.i_vert_2.AND.NT(i_tetra1,4).NE.i_vert_3)  i_vert_opp = NT (i_tetra1,4)
               GOTO 100
           END IF
       END IF
   END IF
END DO
100 CONTINUE
DO i=1,NTETRA
   IF (i_vert_1.EQ.NT(i,1).OR.i_vert_1.EQ.NT(i,2).OR.i_vert_1.EQ.NT(i,3).OR.i_vert_1.EQ.NT(i,4)) THEN
       IF (i_vert_2.EQ.NT(i,1).OR.i_vert_2.EQ.NT(i,2).OR.i_vert_2.EQ.NT(i,3).OR.i_vert_2.EQ.NT(i,4)) THEN
           IF (i_vert_3.EQ.NT(i,1).OR.i_vert_3.EQ.NT(i,2).OR.i_vert_3.EQ.NT(i,3).OR.i_vert_3.EQ.NT(i,4)) THEN
               TETRA_IN_STRUCTURE (i) = .false.
               IF (i.ne.i_tetra1) i_tetra2 = i
           END IF
       END IF
   END IF
END DO
END SUBROUTINE 
 
 
 
!*************************************************************************
!*************************************************************************
 
 
 
SUBROUTINE  CONTROLLA_NODO_PRIMARIO_PIU_VICINO (xs,ys,zs,distanza_minima,num_nodo)
USE CONFIG_MESH
IMPLICIT NONE
DOUBLE PRECISION :: xs,ys,zs
DOUBLE PRECISION :: distanza_minima
DOUBLE PRECISION :: lunghezza_lato
INTEGER :: i
INTEGER :: num_nodo
num_nodo = 1
distanza_minima = DSQRT(((X(1)-Xs)**2.0 + (Y(1)-Ys)**2.0 + (Z(1)-Zs)**2.0))
DO i=1 , NNOD_PRIMARI
    lunghezza_lato = DSQRT(((X(i)-Xs)**2.0 + (Y(i)-Ys)**2.0 + (Z(i)-Zs)**2.0))
    IF (lunghezza_lato.le.distanza_minima) THEN
        distanza_minima = lunghezza_lato
        num_nodo = i
    END IF
END DO
END SUBROUTINE
 
 
 
!*************************************************************************
!*************************************************************************
 
 
 
SUBROUTINE SPOSTA_NODI_SECONDARI_SUI_NODI_PRIMARI
USE CONFIG_MESH
USE CONFIG_PARAMETRI
IMPLICIT NONE
DOUBLE PRECISION :: lato_minimo_presente
DOUBLE PRECISION :: distanza
INTEGER          :: num_nodo
INTEGER          :: in
lato_minimo_presente = DMax1 (XMAX-XMIN,YMAX-YMIN,ZMAX-ZMIN) / FATTORE_REFINIMENT
WRITE (1000,*) 'Lato_Minimo_Presente',lato_minimo_presente
WRITE (1000,*) 'NNOD_PRIMARI ',NNOD_PRIMARI
 
DO In = NNOD_PRIMARI + 1  ,NNOD
    CALL CONTROLLA_NODO_PRIMARIO_PIU_VICINO (x(IN),y(IN),z(IN),distanza,num_nodo)
    WRITE (1000,*) 'distanza',distanza
    IF (distanza.LT.(lato_minimo_presente/15.0)) THEN
         WRITE (1000,*) 'Scambio Effettuato'
         X (in) = x (num_nodo)
         Y (in) = y (num_nodo)
         Z (in) = z (num_nodo)
    END IF
END DO
END SUBROUTINE
 
 
 
!*************************************************************************
!*************************************************************************
 
 
 
SUBROUTINE MIGLIORA_QUALITA_MUOVENDO_LATI_SECONDARI
USE CONFIG_MESH
USE CONFIG_PARAMETRI
USE WINTERACTER
IMPLICIT NONE
TYPE (win_message) :: MESSAGE
INTEGER :: itype
INTEGER :: I,J,k
DOUBLE PRECISION :: x_aux,y_aux,z_aux
DOUBLE PRECISION :: qualita
DOUBLE PRECISION :: volume_minimo_per_cento,volume_parallalelepipedo
DOUBLE PRECISION :: volume_minimo_old
INTEGER          :: i_minimo
CALL TROVA_TETRAEDRO_VOLUME_MINIMO(i_minimo,volume_minimo_per_cento,volume_parallalelepipedo)
volume_minimo_old = volume_minimo_per_cento
WRITE (1000,*) 'Nodi primari', NNOD_PRIMARI
WRITE (1000,*) 'Nodi        ', NNOD
DO i = NNOD_PRIMARI + 1, NNOD
    CALL STAMPA_INFORMAZIONI_SU_SCHERMO(0,0)
    CALL CALCOLO_COEFFICIENTE_QUALITA_GLOBALE ()
    QUALITA = QUALITA_MEDIA
    IF (X(i).eq.XMIN.OR.X(i).eq.XMAX) GOTO 2000
    IF (Y(i).eq.YMIN.OR.Y(i).eq.YMAX) GOTO 2000
    IF (Z(i).eq.ZMIN.OR.Z(i).eq.ZMAX) GOTO 2000
    x_aux = x(i)
    y_aux = y(i)
    z_aux = z(i)
    DO j=1,NTETRA
         CALL WMESSAGEPEEK (ITYPE, MESSAGE)
         IF (ITYPE.eq.CLOSEREQUEST) STOP
         IF (i.EQ.NT(J,1).OR.i.EQ.NT(J,2).OR.i.EQ.NT(J,3).OR.i.EQ.NT(J,4)) THEN
              DO k = 1,4
                   IF (i.NE.NT(j,k)) THEN
                         WRITE (1000,*) 'Numero nodo  ',i
                         X(i) = (X_aux * 8 + X(NT(j,k))) / 9.0
                         Y(i) = (Y_aux * 8 + Y(NT(j,k))) / 9.0
                         Z(i) = (Z_aux * 8 + Z(NT(j,k))) / 9.0
                         WRITE (1000,*) 'x_aux, y_aux, z_aux  ',x_aux,y_aux,z_aux
                         WRITE (1000,*) 'nt                   ',X(NT(j,k)),Y(NT(j,k)),Z(NT(j,k))
                         WRITE (1000,*) 'x(i) y(i) z(i)       ',X(i) , Y(i) , Z(i)
                         CALL TROVA_TETRAEDRO_VOLUME_MINIMO(i_minimo,volume_minimo_per_cento,volume_parallalelepipedo)
                         WRITE (1000,*) 'Volume minimo per cento', volume_minimo_per_cento
                         WRITE (1000,*) 'volume minimo old      ', volume_minimo_old
                         IF (volume_minimo_per_cento.lt.volume_minimo_old) THEN
                               x(i) = x_aux
                               y(i) = y_aux
                               z(i) = z_aux
                               WRITE (1000,*) '----- Volume minimo peggiorato -----'
                               GOTO 100
                         ELSE
                               volume_minimo_old = volume_minimo_per_cento
                         END IF
                         CALL CALCOLO_COEFFICIENTE_QUALITA_GLOBALE ()
                         WRITE (1000,*) 'Qualita_media',QUALITA_MEDIA
                         IF (QUALITA_MEDIA.lt.QUALITA) THEN
                               x(i) = x_aux
                               y(i) = y_aux
                               z(i) = z_aux
                         ELSE
                               QUALITA = QUALITA_MEDIA
                         ENDIF
                         100 CONTINUE
                   END IF
              END DO
         END IF
    END DO
2000 CONTINUE
END DO
CALL WDIALOGHIDE
END SUBROUTINE
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE  TROVA_LATO_PIU_LUNGO_CON_REGOLE (it,i,j,lunghezza_best,lunghezza_worse,i_best,j_best,it_best)
USE CONFIG_MESH 
IMPLICIT NONE
DOUBLE PRECISION :: lunghezza_best,lunghezza_lato
DOUBLE PRECISION :: lunghezza_worse
INTEGER :: it,i,j,i_best,j_best,it_best
lunghezza_lato = DSQRT(((X(NT(It,i))-X(NT(It,j)))**2.0 + (Y(NT(It,i))-Y(NT(It,j)))**2.0 + (Z(NT(It,i))-Z(NT(It,j)))**2.0))
IF (lunghezza_lato.gt.lunghezza_best.and.lunghezza_lato.lt.lunghezza_worse-1.E-2)  then
   lunghezza_best=lunghezza_lato
   i_best = i
   j_best = j
   it_best = it
END IF
END SUBROUTINE
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE  TROVA_LATO_PIU_LUNGO (it,i,j,lunghezza_best,i_best,j_best,it_best)
USE CONFIG_MESH 
IMPLICIT NONE
DOUBLE PRECISION ::   lunghezza_best,lunghezza_lato
INTEGER :: it,i,j,i_best,j_best,it_best
lunghezza_lato = DSQRT(((X(NT(It,i))-X(NT(It,j)))**2.0 + (Y(NT(It,i))-Y(NT(It,j)))**2.0 + (Z(NT(It,i))-Z(NT(It,j)))**2.0))
IF (lunghezza_lato.gt.lunghezza_best)  then
   lunghezza_best=lunghezza_lato
   i_best = i
   j_best = j
   it_best = it
END IF
END SUBROUTINE
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE RIFINISCI_LA_MESH_AL_MASSIMO
USE CONFIG_MESH
USE CONFIG_PARAMETRI
USE CONFIG_FEM
USE CONFIG_ZERI
USE WINTERACTER
USE IDENTIFIERS
IMPLICIT NONE
TYPE (win_message) :: MESSAGE
INTEGER :: ITYPE 
LOGICAL :: continua
INTEGER :: ind
INTEGER :: k_div
INTEGER :: It
INTEGER :: k_migliore
INTEGER :: i_best,j_best,it_best
DOUBLE PRECISION :: lunghezza_best
DOUBLE PRECISION :: lato_minimo_parall
DOUBLE PRECISION :: x1,y1,z1,x2,y2,z2
CALL WDIALOGSELECT     (IDD_DIALOG001)
CALL CALCOLO_DEL_FATTORE_DI_SCALA_PER_LE_PROIEZIONI_ORTOGONALI ()
CALL WMESSAGEBOX (0,0,1,'Questa subroutine deve essere stoppata','')
lato_minimo_parall = DMax1 (XMAX-XMIN,YMAX-YMIN,ZMAX-ZMIN)
WRITE (1000,*) 'lato_minimo_parall',lato_minimo_parall
lunghezza_best = lato_minimo_parall
FATTORE_REFINIMENT = FATTORE_REFINIMENT + 0.5
WRITE (1000,*) 'FATTORE_REFINIMENT',FATTORE_REFINIMENT
continua = .TRUE.
DO WHILE (CONTINUA)
continua =.FALSE.
   DO It=1, NTETRA
!      IF (MATERIALE_TETRAEDRO (it).EQ.1) goto 1000
 
       lunghezza_best = 0.0
       CALL TROVA_LATO_PIU_LUNGO (it,1,2,lunghezza_best,i_best,j_best,it_best)
       CALL TROVA_LATO_PIU_LUNGO (it,1,3,lunghezza_best,i_best,j_best,it_best)
       CALL TROVA_LATO_PIU_LUNGO (it,1,4,lunghezza_best,i_best,j_best,it_best)
       CALL TROVA_LATO_PIU_LUNGO (it,2,3,lunghezza_best,i_best,j_best,it_best)
       CALL TROVA_LATO_PIU_LUNGO (it,2,4,lunghezza_best,i_best,j_best,it_best)
       CALL TROVA_LATO_PIU_LUNGO (it,3,4,lunghezza_best,i_best,j_best,it_best)
!       WRITE (1000,*) 'lunghezza_best',lunghezza_best
       IF (lunghezza_best.GE.(lato_minimo_parall)/FATTORE_REFINIMENT) THEN
           continua = .TRUE.
           x1 = X(NT(it,i_best))
           y1 = Y(NT(it,i_best))
           z1 = Z(NT(it,i_best))
           x2 = X(NT(it,j_best))
           y2 = Y(NT(it,j_best))
           z2 = Z(NT(it,j_best))
           k_div = 80
           k_migliore = k_div / 2
           NNOD = NNOD + 1
           X (NNOD) = x1 + (x2-x1) * k_migliore / k_div
           Y (NNOD) = y1 + (y2-y1) * k_migliore / k_div
           Z (NNOD) = z1 + (z2-z1) * k_migliore / k_div
           CALL POSIZIONE_PUNTO_IE_RISPETTO_TETRAEDRO (it, X(NNOD), y(NNOD), Z(NNOD), ind, ZERO_NODO_IN_TETRAEDRO)
           CALL IMMISSIONE_NODO_NEL_TETRAEDRO (it,NNOD,ind)
           CONDIZIONE_INIZIALE (NNOD) = -1.0
           CALL DISEGNA_LE_PROIEZIONI_ORTOGONALI_X_CARCIOFO (0)
           CALL DISEGNA_LE_PROIEZIONI_ORTOGONALI_DEI_TETRAEDRI ()
           CALL WDIALOGPUTINTEGER (idf_integer3 , NTETRA)
           CALL WDIALOGPUTINTEGER (idf_integer2 , NLAT)
           CALL WDIALOGPUTINTEGER (idf_integer1 , NNOD)
           CALL WDIALOGPUTREAL    (IDF_REAL8, REAL (FATTORE_REFINIMENT))
           CALL WMESSAGEPEEK      (ITYPE, MESSAGE)
           IF (ITYPE.eq.CLOSEREQUEST) STOP
       END IF
       1000 CONTINUE
   END DO
END DO
PAUSE
CALL WDIALOGHIDE
END SUBROUTINE
 
 
 
!*************************************************************************
!*************************************************************************
 
 
 
SUBROUTINE CALCOLA_QUALITA_PEGGIORE_DEL_TETRAEDRO (qualita_peggiore)
USE CONFIG_MESH
IMPLICIT NONE
INTEGER :: I
INTEGER :: It
INTEGER :: i_best,j_best,k_best
INTEGER :: i_tetra_peggiore
DOUBLE PRECISION :: x1,y1,z1,x2,y2,z2,x3,y3,z3
DOUBLE PRECISION :: area_best
DOUBLE PRECISION :: qualita_peggiore
qualita_peggiore = 1.0
DO It = 1, NTETRA
    CALL TROVA_TETRA_CON_QUALITA_PEGGIORE (it,qualita_peggiore,i_tetra_peggiore)
END DO
WRITE (1000,*) 'qualita_peggiore',i_tetra_peggiore,qualita_peggiore
END SUBROUTINE
 
 
 
!*************************************************************************
!*************************************************************************
 
 
 
SUBROUTINE  TROVA_TETRA_CON_QUALITA_PEGGIORE (i,qualita_peggiore,itetra_peggiore)
USE CONFIG_MESH
USE CONFIG_PARAMETRI
IMPLICIT NONE
INTEGER :: i
DOUBLE PRECISION :: qt
DOUBLE PRECISION :: qualita_tetra
DOUBLE PRECISION :: qualita_peggiore
INTEGER :: itetra_peggiore
qualita_tetra = QT (I)
IF (qualita_tetra.LE.qualita_peggiore) THEN
     qualita_peggiore = qualita_tetra
     itetra_peggiore = I
END IF
END SUBROUTINE
 
 
 
!*************************************************************************
!*************************************************************************
 
SUBROUTINE RIFINISCI_LA_MESH_AL_MASSIMO_IN_BASE_ALLA_QUALITA
USE CONFIG_MESH
USE CONFIG_PARAMETRI
IMPLICIT NONE
INTEGER :: I,k,k_div
INTEGER :: It
INTEGER :: k_migliore
INTEGER :: i_best,j_best
INTEGER :: i_tetra_peggiore
DOUBLE PRECISION :: lunghezza_best
DOUBLE PRECISION :: x1,y1,z1,x2,y2,z2
DOUBLE PRECISION :: xs,ys,zs
DOUBLE PRECISION :: distanza_migliore,distanza_minima
DOUBLE PRECISION :: qualita_peggiore
lunghezza_best = yMAX-yMIN
FATTORE_REFINIMENT = FATTORE_REFINIMENT + 0.2
DO I= 1, 10
   qualita_peggiore = 1.0
   DO It = 1, NTETRA
        CALL TROVA_TETRA_CON_QUALITA_PEGGIORE (it,qualita_peggiore,i_tetra_peggiore)
   END DO
   WRITE (1000,*) 'qualita_peggiore',i_tetra_peggiore,qualita_peggiore
   lunghezza_best = 0.0
   CALL TROVA_LATO_PIU_LUNGO (i_tetra_peggiore,1,2,lunghezza_best,i_best,j_best,i_tetra_peggiore)
   CALL TROVA_LATO_PIU_LUNGO (i_tetra_peggiore,1,3,lunghezza_best,i_best,j_best,i_tetra_peggiore)
   CALL TROVA_LATO_PIU_LUNGO (i_tetra_peggiore,1,4,lunghezza_best,i_best,j_best,i_tetra_peggiore)
   CALL TROVA_LATO_PIU_LUNGO (i_tetra_peggiore,2,3,lunghezza_best,i_best,j_best,i_tetra_peggiore)
   CALL TROVA_LATO_PIU_LUNGO (i_tetra_peggiore,2,4,lunghezza_best,i_best,j_best,i_tetra_peggiore)
   CALL TROVA_LATO_PIU_LUNGO (i_tetra_peggiore,3,4,lunghezza_best,i_best,j_best,i_tetra_peggiore)
   x1 = X(NT(i_tetra_peggiore,i_best))
   y1 = Y(NT(i_tetra_peggiore,i_best))
   z1 = Z(NT(i_tetra_peggiore,i_best))
   x2 = X(NT(i_tetra_peggiore,j_best))
   y2 = Y(NT(i_tetra_peggiore,j_best))
   z2 = Z(NT(i_tetra_peggiore,j_best))
   distanza_migliore = 0.0
   k_div = 40
   DO k= 18, k_div-18
       xs = x1 + (x2-x1) * k / k_div
       ys = y1 + (y2-y1) * k / k_div
       zs = z1 + (z2-z1) * k / k_div
       CALL CONTROLLO_DISTANZA_MINIMA_CON_ALTRI_NODI (xs,ys,zs,distanza_minima)
       IF (distanza_minima.gt.distanza_migliore) THEN
             distanza_migliore = distanza_minima
             k_migliore = k
       END IF 
   END DO
   NNOD = NNOD + 1
   k_migliore=20
   X (NNOD) = x1 + (x2-x1) * k_migliore / k_div
   Y (NNOD) = y1 + (y2-y1) * k_migliore / k_div
   Z (NNOD) = z1 + (z2-z1) * k_migliore / k_div
   CALL TETRAEDIZZAZIONE_PER_NODI (NNOD, 0, 0, NNOD)
END DO 
CALL WDIALOGHIDE
END SUBROUTINE
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE RIFINISCI_LA_MESH_AL_MASSIMO_IN_BASE_ALLA_QUALITA_CON_BARICENTRO
USE CONFIG_MESH
USE CONFIG_PARAMETRI
IMPLICIT NONE
INTEGER :: I
INTEGER :: It
INTEGER :: i_tetra_peggiore
DOUBLE PRECISION :: x1,y1,z1,x2,y2,z2,x3,y3,z3,x4,y4,z4
DOUBLE PRECISION :: qualita_peggiore
FATTORE_REFINIMENT = FATTORE_REFINIMENT + 0.2
DO I= 1, 10
   qualita_peggiore = 1.0
   DO It = 1, NTETRA
        CALL TROVA_TETRA_CON_QUALITA_PEGGIORE (it,qualita_peggiore,i_tetra_peggiore)
   END DO
   WRITE (1000,*) 'qualita_peggiore',i_tetra_peggiore,qualita_peggiore
   x1 = X(NT(i_tetra_peggiore,1))
   y1 = Y(NT(i_tetra_peggiore,1))
   z1 = Z(NT(i_tetra_peggiore,1))
   x2 = X(NT(i_tetra_peggiore,2))
   y2 = Y(NT(i_tetra_peggiore,2))
   z2 = Z(NT(i_tetra_peggiore,2))
   x3 = X(NT(i_tetra_peggiore,3))
   y3 = Y(NT(i_tetra_peggiore,3))
   z3 = Z(NT(i_tetra_peggiore,3))
   x4 = X(NT(i_tetra_peggiore,4))
   y4 = Y(NT(i_tetra_peggiore,4))
   z4 = Z(NT(i_tetra_peggiore,4))
   NNOD = NNOD + 1
   X (NNOD) = (x1 + x2 + x3 + x4) / 4.0
   Y (NNOD) = (y1 + y2 + y3 + y4) / 4.0
   Z (NNOD) = (z1 + z2 + z3 + z4) / 4.0
   CALL TETRAEDIZZAZIONE_PER_NODI (NNOD, 0, 0, NNOD)
END DO 
CALL WDIALOGHIDE
END SUBROUTINE
 
 
!*************************************************************************
!*************************************************************************
 
SUBROUTINE RIFINISCI_LA_MESH_AL_MASSIMO_IN_BASE_ALLA_QUALITA_CON_AREA
USE CONFIG_MESH
IMPLICIT NONE
INTEGER :: I
INTEGER :: It
INTEGER :: i_best,j_best,k_best
INTEGER :: i_tetra_peggiore
DOUBLE PRECISION :: x1,y1,z1,x2,y2,z2,x3,y3,z3
DOUBLE PRECISION :: area_best
DOUBLE PRECISION :: qualita_peggiore
DO I= 1, 10
   qualita_peggiore = 1.0
   DO It = 1, NTETRA
        CALL TROVA_TETRA_CON_QUALITA_PEGGIORE (it,qualita_peggiore,i_tetra_peggiore)
   END DO
   WRITE (1000,*) 'qualita_peggiore',i_tetra_peggiore,qualita_peggiore
   area_best = 0.0
   CALL TROVA_AREA_PIU_GRANDE (i_tetra_peggiore,1,2,3,area_best,i_best,j_best,k_best)
   CALL TROVA_AREA_PIU_GRANDE (i_tetra_peggiore,1,2,4,area_best,i_best,j_best,k_best)
   CALL TROVA_AREA_PIU_GRANDE (i_tetra_peggiore,2,3,4,area_best,i_best,j_best,k_best)
   CALL TROVA_AREA_PIU_GRANDE (i_tetra_peggiore,1,3,4,area_best,i_best,j_best,k_best)
   x1 = X(NT(i_tetra_peggiore,i_best))
   y1 = Y(NT(i_tetra_peggiore,i_best))
   z1 = Z(NT(i_tetra_peggiore,i_best))
   x2 = X(NT(i_tetra_peggiore,j_best))
   y2 = Y(NT(i_tetra_peggiore,j_best))
   z2 = Z(NT(i_tetra_peggiore,j_best))
   x3 = X(NT(i_tetra_peggiore,k_best))
   y3 = Y(NT(i_tetra_peggiore,k_best))
   z3 = Z(NT(i_tetra_peggiore,k_best))
   NNOD = NNOD + 1
   X (NNOD) = (x1 + x2 + x3) / 3.0
   Y (NNOD) = (y1 + y2 + y3) / 3.0
   Z (NNOD) = (z1 + z2 + z3) / 3.0
   CALL TETRAEDIZZAZIONE_PER_NODI (NNOD, 0, 0, NNOD)
END DO 
CALL WDIALOGHIDE
END SUBROUTINE
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE TROVA_AREA_PIU_GRANDE (i_tetra,i,j,k,area_best,i_best,j_best,k_best)
USE CONFIG_MESH
IMPLICIT NONE
INTEGER :: i,j,k
INTEGER :: i_best,j_best,k_best
INTEGER :: i_tetra
DOUBLE PRECISION :: area_best,area
CALL CALCOLO_AREA_DEL_TRIANGOLO (NT(i_tetra,i), NT(i_tetra,j), NT(i_tetra,k), area)
IF (area.gt.area_best) THEN
    i_best = i
    j_best = j
    k_best = k
    area_best = area
END IF 
END SUBROUTINE 
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE CREA_NODI_SECONDARI_PER_MIGLIORARE
USE CONFIG_MESH
USE CONFIG_PARAMETRI
IMPLICIT NONE
INTEGER :: i,k,k_div
INTEGER :: k_migliore
INTEGER :: It
INTEGER :: i_best,j_best,it_best
DOUBLE PRECISION :: lunghezza_best
DOUBLE PRECISION :: x1, y1, z1, x2, y2, z2
DOUBLE PRECISION :: xs, ys, zs
DOUBLE PRECISION :: distanza_migliore
DOUBLE PRECISION :: distanza_minima
DO I=1,10
   lunghezza_best = 0.0
   DO It=1, NTETRA
       CALL TROVA_LATO_PIU_LUNGO (it,1,2,lunghezza_best,i_best,j_best,it_best)
       CALL TROVA_LATO_PIU_LUNGO (it,1,3,lunghezza_best,i_best,j_best,it_best)
       CALL TROVA_LATO_PIU_LUNGO (it,1,4,lunghezza_best,i_best,j_best,it_best)
       CALL TROVA_LATO_PIU_LUNGO (it,2,3,lunghezza_best,i_best,j_best,it_best)
       CALL TROVA_LATO_PIU_LUNGO (it,2,4,lunghezza_best,i_best,j_best,it_best)
       CALL TROVA_LATO_PIU_LUNGO (it,3,4,lunghezza_best,i_best,j_best,it_best)
   END DO
   x1 = X(NT(it_best,i_best))
   y1 = Y(NT(it_best,i_best))
   z1 = Z(NT(it_best,i_best))
   x2 = X(NT(it_best,j_best))
   y2 = Y(NT(it_best,j_best))
   z2 = Z(NT(it_best,j_best))
   distanza_migliore = 0.0
   k_div = 40
   DO k= 14, k_div-14
       xs = x1 + (x2-x1) * k / k_div
       ys = y1 + (y2-y1) * k / k_div
       zs = z1 + (z2-z1) * k / k_div
       CALL CONTROLLO_DISTANZA_MINIMA_CON_ALTRI_NODI (xs,ys,zs,distanza_minima)
       IF (distanza_minima.gt.distanza_migliore) THEN
             distanza_migliore = distanza_minima
             k_migliore = k
       END IF 
   END DO
   NNOD = NNOD + 1
   X (NNOD) = x1 + (x2-x1) * k_migliore / k_div
   Y (NNOD) = y1 + (y2-y1) * k_migliore / k_div
   Z (NNOD) = z1 + (z2-z1) * k_migliore / k_div
   CALL TETRAEDIZZAZIONE_PER_NODI (NNOD, 1, 0, NNOD)
END DO
CALL WDIALOGHIDE
END SUBROUTINE
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE  CONTROLLO_DISTANZA_MINIMA_CON_ALTRI_NODI (xs,ys,zs,distanza_minima)
USE CONFIG_MESH
IMPLICIT NONE
DOUBLE PRECISION :: xs,ys,zs
DOUBLE PRECISION :: distanza_minima
DOUBLE PRECISION :: lunghezza_lato
INTEGER :: i
distanza_minima = DSQRT(((X(1)-Xs)**2.0 + (Y(1)-Ys)**2.0 + (Z(1)-Zs)**2.0))
DO i=1,NNOD
    lunghezza_lato = DSQRT(((X(i)-Xs)**2.0 + (Y(i)-Ys)**2.0 + (Z(i)-Zs)**2.0))
    IF (lunghezza_lato.le.distanza_minima) distanza_minima = lunghezza_lato
END DO
DO i = 1, 8
    lunghezza_lato = DSQRT(((X(MaxNODI + i)-Xs)**2.0 + (Y(MaxNODI + i)-Ys)**2.0 + (Z(MaxNODI + i)-Zs)**2.0))
    IF (lunghezza_lato/1.0.le.distanza_minima) distanza_minima = lunghezza_lato
END DO
END SUBROUTINE
 
 
 
!*************************************************************************
!*************************************************************************
 
 
SUBROUTINE RIFINISCI_LA_MESH_CON_NODI_SECONDARI_LONTANI
USE CONFIG_MESH
USE CONFIG_PARAMETRI
 
IMPLICIT NONE
LOGICAL :: continua
INTEGER :: k_div
INTEGER :: It,K
INTEGER :: k_migliore
INTEGER :: i_best,j_best,it_best
DOUBLE PRECISION :: lunghezza_best
DOUBLE PRECISION :: lato_minimo_parall
DOUBLE PRECISION :: x1,y1,z1,x2,y2,z2
DOUBLE PRECISION :: distanza_migliore,distanza_minima
DOUBLE PRECISION :: xs,ys,zs
 
lato_minimo_parall = DMax1 (XMAX-XMIN,YMAX-YMIN,ZMAX-ZMIN)
WRITE (1000,*) 'lato_minimo_parall',lato_minimo_parall
lunghezza_best = lato_minimo_parall
FATTORE_REFINIMENT = FATTORE_REFINIMENT + 0.3
WRITE (1000,*) 'FATTORE_REFINIMENT',FATTORE_REFINIMENT
continua = .TRUE.
DO WHILE (CONTINUA)
continua =.FALSE.
   DO It=1, NTETRA
       lunghezza_best = 0.0
       CALL TROVA_LATO_PIU_LUNGO (it,1,2,lunghezza_best,i_best,j_best,it_best)
       CALL TROVA_LATO_PIU_LUNGO (it,1,3,lunghezza_best,i_best,j_best,it_best)
       CALL TROVA_LATO_PIU_LUNGO (it,1,4,lunghezza_best,i_best,j_best,it_best)
       CALL TROVA_LATO_PIU_LUNGO (it,2,3,lunghezza_best,i_best,j_best,it_best)
       CALL TROVA_LATO_PIU_LUNGO (it,2,4,lunghezza_best,i_best,j_best,it_best)
       CALL TROVA_LATO_PIU_LUNGO (it,3,4,lunghezza_best,i_best,j_best,it_best)
       IF (lunghezza_best.GE.(lato_minimo_parall)/FATTORE_REFINIMENT) THEN
           continua = .TRUE.
           x1 = X(NT(it,i_best))
           y1 = Y(NT(it,i_best))
           z1 = Z(NT(it,i_best))
           x2 = X(NT(it,j_best))
           y2 = Y(NT(it,j_best))
           z2 = Z(NT(it,j_best))
           k_div = 80
           k_migliore = k_div / 2
           DO k= 35, k_div-35
               distanza_migliore = 0.0
               xs = x1 + (x2-x1) * k / k_div
               ys = y1 + (y2-y1) * k / k_div
               zs = z1 + (z2-z1) * k / k_div
               CALL CONTROLLO_DISTANZA_MINIMA_CON_ALTRI_NODI (xs,ys,zs,distanza_minima)
               WRITE (1000,*) 'distanza_minima ',distanza_minima
               IF (distanza_minima.gt.distanza_migliore) THEN
                   distanza_migliore = distanza_minima
                   k_migliore = k
               END IF
           END DO
           NNOD = NNOD + 1
           X (NNOD) = x1 + (x2-x1) * k_migliore / k_div
           Y (NNOD) = y1 + (y2-y1) * k_migliore / k_div
           Z (NNOD) = z1 + (z2-z1) * k_migliore / k_div
           CALL TETRAEDIZZAZIONE_PER_NODI (NNOD, 0, 0, NNOD)
       END IF
   END DO
END DO
CALL WDIALOGHIDE
END SUBROUTINE

Nelle Figure 4-25, 4-26 viene mostrato un modello di un polo di una macchina sincrona a flusso radiale con 3 cave per polo per fase. Il modello prevede una realizzazione schematica delle testate di macchina L’eccitazione viene prodotta da un magnete permanente montato sulla parte rotorica della macchina. Il processo di reticolazione segue l’evoluzione mostrata nelle Figure 4-27, 28, 29, 30, 31, 32, 33, 34. Un particolare della mesh dopo il processo di refinement viene riportata per gli avvolgimenti di macchina nelle Figure 4-35, 36.

Figura 425   Struttura globale della macchina radiale. Viene modellato un polo elettrico di una macchina ad 8 poli

Figura – Particolare della struttura degli avvolgimenti di macchina

       

Figura 427   Wireframe e superficie della mesh, dopo l’inserimento di 500 lati

Figura 428   Wireframe e superficie della mesh, dopo l’inserimento di 1000 lati

Figura 429   Wireframe e superficie della mesh, dopo l’inserimento di 1500 lati

Figura 430   Wireframe e superficie della mesh, dopo l’inserimento di 2000 lati

Figura 431   Wireframe e superficie della mesh, dopo l’inserimento di 2500 lati

Figura 432   Wireframe e superficie della mesh, dopo l’inserimento di 3000 lati

Figura 433   Wireframe e superficie della mesh, dopo l’inserimento di 3500 lati

Figura 434   Risultato finale del processo di refinement

Figura 435   Vista di due avvolgimenti statorici con particolare della discretizzazione spaziale prodotto dal processo di triangolazione

Figura 436   Vista di un avvolgimento statorico di fase. La mesh è stata già rifinita con il criterio frattale

 

  1. BIBLIOGRAFIA

  1. Farin, G. – Curves and Surfaces for Computer Aided Geometric Design -. Academic Press Inc. Ltd: London, I Edition, 1988.

  2. Borouchaki, H., Hecht, F., Saltel, E. & George, P.L. – Reasonably Efficient Delaunay Based Mesh Generator in 3 Dimensions – . Proc. of 4th International Meshing Roundtable, Sandia National Laboratories, pp. 3-14, 1995.

  3. Freitag, L.A. – Tetrahedral Mesh Improvement Using Swapping and Smoothing – International Journal for Numerical Methods in Engineering, Wiley, 40, pp. 3979-4002, 1997.

  4. Freitag, L.A. & Ollivier­Gooch, C. – A Comparison of Tetrahedral Mesh Improvement Techniques – Proceedings of 5th International Meshing Roundtable, Sandia National Laboratories, pp. 87-100, 1996.

  5. George, P.L., Hecht, F. & Saltel, E. – Automatic mesh generator with specified boundary – Comput. Methods Appl. Mech. Engrg, 92, pp. 269-288, 1991.

  6. Hazlewood, C. – Approximating constrained tetrahedralizations – Computer Aided Geometric Design, 10, pp. 67-87, 1993.

  7. Golias, N.A. & Tsiboukis, T.D. – An Approach to Refining Three­Dimensional Tetrahedral Meshes Based on Delaunay Transformations – Int. J. for Numerical Methods in Engineering, 37, pp. 793-812, 1994.

  8. Borouchaki, H., Hecht, F., Saltel, E. & George, P.L. – Reasonably Efficient Delaunay Based Mesh Generator in 3 Dimensions – Proc. of 4th International Meshing Roundtable, Sandia National Laboratories, pp. 3-14, 1995.

  9. Yamashita, H.; Yamaji, A.; Cingoski, V.; Kaneda, K.;Magnetics – A novel tetrahedral mesh generation method for rotating machines including end-coil region – IEEE Transactions on , Volume: 32 , Issue: 3 , May 1996 Pages:1353 – 1356.

  10. Cendes, Z.; Shenton, D. – Adaptive mesh refinement in the finite element computation of magnetic fields – Magnetics, IEEE Transactions on , Volume: 21 , Issue: 5 , Sep 1985 Pages:1811 – 1816.

  11. Kahler, K.; Haber, J.; Seidel, H.-P – Dynamic refinement of deformable triangle meshes for rendering – Computer Graphics International 2001. Proceedings , 3-6 July 2001 Pages:285 – 290.

  12. Ivrissimtzis, L.; Seidel, H.-P – Polyhedra operators for mesh refinement -Geometric Modeling and Processing, 2002. Proceedings , 10-12 July 2002 Pages:132 – 137.

  13. O’Dwyer, J.; Evans, P. – Triangular element refinement in automatic adaptive mesh generation – Magnetics, IEEE Transactions on , Volume: 33 , Issue: 2, March 1997 Pages:1740 – 1743.

  14. Adaptive refinement strategies in three dimensions Golias, N.A.; Tsiboukis, T.D.;Magnetics, IEEE Transactions on , Volume: 29 , Issue: 2 , Mar 1993  Pages:1886 – 1889.

  15. Real-time refinement and simplification of adaptive triangular meshesVolkov, V.; Ling Li; Visualization, 2003. VIS 2003. IEEE , 19-24 Oct. 2003
    Pages:155 – 162.

  16.  Selective refinement queries for volume visualization of unstructured tetrahedral meshesCignoni, P.; De Floriani, L.; Magillo, P.; Puppo, E.; Scopigno, R.;Visualization and Computer Graphics, IEEE Transactions on , Volume: 10 , Issue: 1 , Jan-Feb 2004 Pages:29 – 45.

  17. Olszewski, P.; Nakata, T.; Takahashi, N.; ujiwara, K. – A simple algorithm for adaptive refinement of tetrahedral meshes combined with edge elements – Magnetics, IEEE Transactions on , Volume: 29 , Issue: 2 , Mar 1993 Pages:1898 – 1901.

  18. Raizer, A., Meunier, G., Coulomb, J.-L – An approach for automatic adaptive mesh refinement in finite element computation of magnetic fields – Magnetics, IEEE Transactions on , Volume: 25 , Issue: 4 , July 1989 Pages:2965 – 2967.

  19. D. Bai, A. Brandt – Local mesh refinement multilevel techniques – SIAM J. of Stat. Comp., 1987.

  20. R.E. Bank, A.H. Sherman – Refinement algorithms and data structure for regular local mesh refinement – IMACS, North Holland, 1983.

  21. Harris, D. – An exponentiation unit for an OpenGL lighting engine – Computers, IEEE Transactions on , Volume: 53 , Issue: 3 , March 2004 Pages:251 – 258.

Bhanirantka, P.; Demange, Y. – OpenGL Volumizer: a toolkit for high quality volume rendering of large data sets – Volume Visualization and Graphics, 2002. Proceedings. IEEE / ACM SIGGRAPH Symposium on , 28-29 Oct. 2002 Pages:45 – 53

Lascia un commento

Il tuo indirizzo email non sarà pubblicato. I campi obbligatori sono contrassegnati *