-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathmain.f90
3117 lines (2530 loc) · 110 KB
/
main.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
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
PROGRAM FDS
! Fire Dynamics Simulator, Main Program, Multiple CPU version.
USE PRECISION_PARAMETERS
USE MESH_VARIABLES
USE GLOBAL_CONSTANTS
USE TRAN
USE DUMP
USE READ_INPUT
USE INIT
USE DIVG
USE PRES
USE MASS
USE PART
USE VEGE
USE VELO
USE RAD
USE OUTPUT_DATA
USE MEMORY_FUNCTIONS
USE HVAC_ROUTINES
USE COMP_FUNCTIONS, ONLY : SECOND, WALL_CLOCK_TIME, SHUTDOWN
USE MATH_FUNCTIONS, ONLY : GAUSSJ
USE DEVICE_VARIABLES
USE WALL_ROUTINES
USE FIRE
USE RADCONS, ONLY: TIME_STEP_INCREMENT
USE CONTROL_FUNCTIONS
USE EVAC
USE TURBULENCE, ONLY: NS_ANALYTICAL_SOLUTION,INIT_TURB_ARRAYS,COMPRESSION_WAVE,STRATIFIED_MIXING_LAYER, &
SYNTHETIC_TURBULENCE,SYNTHETIC_EDDY_SETUP,GET_REV_turb
USE COMPLEX_GEOMETRY, ONLY: INIT_IBM,GET_REV_geom,GET_CUTCELL_AREA
USE OPENMP
USE MPI
USE SCRC, ONLY: SCARC_SETUP, SCARC_SOLVER, SCARC_TIMINGS, GET_REV_SCRC
IMPLICIT NONE
! Miscellaneous declarations
CHARACTER(255), PARAMETER :: mainid='$Id: main.f90 $'
CHARACTER(255), PARAMETER :: mainrev='$WFDS based on FDS Revision: 9906 $'
CHARACTER(255), PARAMETER :: maindate='$Date: $'
LOGICAL :: EX,EX_GLOBAL,DIAGNOSTICS,EXCHANGE_RADIATION,EXCHANGE_EVACUATION=.FALSE.,FIRST_PASS
INTEGER :: LO10,NM,IZERO,REVISION_NUMBER,IOS
CHARACTER(255) :: REVISION_DATE
REAL(EB) :: T_MAX,T_MIN
REAL(EB), ALLOCATABLE, DIMENSION(:) :: T,TC_GLB,TC_LOC,DT_SYNC,DT_NEXT_SYNC,TI_LOC,TI_GLB, &
DSUM_ALL,PSUM_ALL,USUM_ALL,DSUM_ALL_LOCAL,PSUM_ALL_LOCAL,USUM_ALL_LOCAL
LOGICAL, ALLOCATABLE, DIMENSION(:,:) :: CONNECTED_ZONES_GLOBAL,CONNECTED_ZONES_LOCAL
INTEGER, ALLOCATABLE, DIMENSION(:) :: MESH_STOP_STATUS
LOGICAL, ALLOCATABLE, DIMENSION(:) :: ACTIVE_MESH,STATE_GLB,STATE_LOC
INTEGER :: NOM,IWW,IW,STOP_STATUS=0
INTEGER, PARAMETER :: N_DROP_ADOPT_MAX=10000
TYPE (MESH_TYPE), POINTER :: M,M4
TYPE (OMESH_TYPE), POINTER :: M2,M3,M5
! MPI stuff
INTEGER :: N,I,IERR=0,STATUS(MPI_STATUS_SIZE)
INTEGER :: BUFFER_SIZE,TAG,PNAMELEN=0,DISP,TAG_EVAC
INTEGER, ALLOCATABLE, DIMENSION(:) :: REQ,PREQ,COUNTS,DISPLS,COUNTS2D,DISPLS2D,COUNTS_TIMERS,DISPLS_TIMERS, &
COUNTS_MASS,DISPLS_MASS,COUNTS_HVAC,DISPLS_HVAC
INTEGER, ALLOCATABLE, DIMENSION(:) :: COUNTS_TREE,DISPLS_TREE
INTEGER :: N_REQ,N_PREQ,N_COMMUNICATIONS
CHARACTER(MPI_MAX_PROCESSOR_NAME) :: PNAME
REAL(EB) :: MPI_WALL_TIME,MPI_WALL_TIME_START
INTEGER, ALLOCATABLE, DIMENSION(:) :: INTEGER_BUFFER_1
REAL(EB), ALLOCATABLE, DIMENSION(:) :: REAL_BUFFER_1
REAL(EB), ALLOCATABLE, DIMENSION(:,:) :: REAL_BUFFER_2,REAL_BUFFER_3,REAL_BUFFER_5,REAL_BUFFER_6,REAL_BUFFER_10
REAL(EB), ALLOCATABLE, DIMENSION(:,:,:) :: REAL_BUFFER_4,REAL_BUFFER_7,REAL_BUFFER_8
REAL(EB), ALLOCATABLE, DIMENSION(:,:,:,:) :: REAL_BUFFER_9
LOGICAL, ALLOCATABLE, DIMENSION(:) :: LOGICAL_BUFFER_1
! Initialize MPI (First executable lines of code)
CALL MPI_INIT(IERR)
CALL MPI_COMM_RANK(MPI_COMM_WORLD, MYID, IERR)
CALL MPI_COMM_SIZE(MPI_COMM_WORLD, NUMPROCS, IERR)
CALL MPI_GET_PROCESSOR_NAME(PNAME, PNAMELEN, IERR)
MPI_WALL_TIME_START = MPI_WTIME()
IF (PNAME/='null') THEN
USE_MPI = .TRUE.
WRITE(LU_ERR,'(A,I3,A,I3,A,A)') 'Process ',MYID,' of ',NUMPROCS-1,' is running on ',PNAME(1:PNAMELEN)
ENDIF
CALL MPI_BARRIER(MPI_COMM_WORLD, IERR)
! Check for OpenMP
CALL OPENMP_CHECK
IF (USE_OPENMP) THEN
IF ( USE_MPI) WRITE(LU_ERR,'(A,A,I3,A)') PNAME(1:PNAMELEN),' has ',OPENMP_AVAILABLE_THREADS,' threads available'
IF (.NOT.USE_MPI) WRITE(LU_ERR,'(I3,A)') OPENMP_AVAILABLE_THREADS,' threads available'
ENDIF
CALL MPI_BARRIER(MPI_COMM_WORLD, IERR)
! Start wall clock timing
WALL_CLOCK_START = WALL_CLOCK_TIME()
! Assign a compilation date (All Nodes)
WRITE(VERSION_STRING,'(A)') '6.0.0'
IF (INDEX(mainrev,':',BACK=.TRUE.)>0) THEN
WRITE(REVISION_DATE,'(A)',IOSTAT=IOS,ERR=5) mainrev(INDEX(mainrev,':')+1:LEN_TRIM(mainrev)-2)
5 REVISION_NUMBER = 0
IF (IOS==0) READ(REVISION_DATE,'(I5)') REVISION_NUMBER
WRITE(REVISION_DATE,'(A)') maindate
CALL GET_REVISION_NUMBER(REVISION_NUMBER,REVISION_DATE)
SVN_REVISION_NUMBER = REVISION_NUMBER
WRITE(COMPILE_DATE,'(A)',IOSTAT=IOS,ERR=10) REVISION_DATE(INDEX(REVISION_DATE,'(')+1:INDEX(REVISION_DATE,')')-1)
10 IF (IOS>0) COMPILE_DATE = 'null'
ENDIF
! Set some constants
CALL SET_OFTEN_USED
! Read input from CHID.data file (All Nodes)
CALL READ_DATA
! Determine if radiation exchanges should be done by default
IF (RADIATION) THEN
EXCHANGE_RADIATION = .TRUE.
ELSE
EXCHANGE_RADIATION = .FALSE.
ENDIF
! Set up send and receive buffer counts and displacements
ALLOCATE(REAL_BUFFER_1(NMESHES))
ALLOCATE(REAL_BUFFER_2(NMESHES,NMESHES))
ALLOCATE(REAL_BUFFER_3(N_TIMERS_DIM,NMESHES))
ALLOCATE(REAL_BUFFER_4(N_TREES_OUT,11,NMESHES))
ALLOCATE(REAL_BUFFER_5(0:MAX_SPECIES,NMESHES))
ALLOCATE(REAL_BUFFER_6(N_DUCTNODES,NMESHES))
ALLOCATE(REAL_BUFFER_7(N_DUCTNODES,N_TRACKED_SPECIES,NMESHES))
ALLOCATE(REAL_BUFFER_8(0:N_ZONE,0:N_ZONE,NMESHES))
ALLOCATE(REAL_BUFFER_9(0:N_ZONE,0:N_ZONE,N_TRACKED_SPECIES,NMESHES))
ALLOCATE(REAL_BUFFER_10(N_DUCTS,NMESHES))
ALLOCATE(INTEGER_BUFFER_1(NMESHES))
ALLOCATE(LOGICAL_BUFFER_1(NMESHES))
ALLOCATE(COUNTS(0:NUMPROCS-1))
ALLOCATE(COUNTS2D(0:NUMPROCS-1))
ALLOCATE(COUNTS_TIMERS(0:NUMPROCS-1))
ALLOCATE(COUNTS_HVAC(0:NUMPROCS-1))
ALLOCATE(COUNTS_MASS(0:NUMPROCS-1))
ALLOCATE(DISPLS(0:NUMPROCS-1))
ALLOCATE(DISPLS2D(0:NUMPROCS-1))
ALLOCATE(DISPLS_MASS(0:NUMPROCS-1))
ALLOCATE(DISPLS_TIMERS(0:NUMPROCS-1))
ALLOCATE(DISPLS_HVAC(0:NUMPROCS-1))
ALLOCATE(COUNTS_TREE(0:NUMPROCS-1))
ALLOCATE(DISPLS_TREE(0:NUMPROCS-1))
COUNTS = 0
DO N=0,NUMPROCS-1
DO NM=1,NMESHES
IF (PROCESS(NM)==N) COUNTS(N) = COUNTS(N) + 1
ENDDO
ENDDO
DISPLS(0) = 0
DO N=1,NUMPROCS-1
DISPLS(N) = COUNTS(N-1) + DISPLS(N-1)
ENDDO
COUNTS2D = COUNTS*NMESHES
DISPLS2D = DISPLS*NMESHES
COUNTS_TIMERS = COUNTS*N_TIMERS_DIM
DISPLS_TIMERS = DISPLS*N_TIMERS_DIM
COUNTS_HVAC = COUNTS*N_DUCTNODES
DISPLS_HVAC = DISPLS*N_DUCTNODES
COUNTS_MASS = COUNTS*(MAX_SPECIES+1)
DISPLS_MASS = DISPLS*(MAX_SPECIES+1)
! Open and write to Smokeview and status file (Master Node Only)
CALL ASSIGN_FILE_NAMES
DO N=0,NUMPROCS-1
IF (MYID==N) CALL WRITE_SMOKEVIEW_FILE
IF (N==NUMPROCS-1) EXIT
CALL MPI_BARRIER(MPI_COMM_WORLD, IERR)
ENDDO
CALL MPI_BARRIER(MPI_COMM_WORLD, IERR)
IF (MYID==0) THEN
OPEN(LU_SMV,FILE=FN_SMV,FORM='FORMATTED', STATUS='OLD',POSITION='APPEND')
CALL WRITE_STATUS_FILES
ENDIF
! Stop all the processes if this is just a set-up run
IF (SET_UP) THEN
CALL INITIALIZE_DIAGNOSTIC_FILE
CALL SHUTDOWN('Stop FDS, Set-up only')
ENDIF
! Set up Time arrays (All Nodes)
ALLOCATE(ACTIVE_MESH(NMESHES),STAT=IZERO)
CALL ChkMemErr('MAIN','ACTIVE_MESH',IZERO)
ACTIVE_MESH = .TRUE.
ALLOCATE(T(NMESHES),STAT=IZERO)
CALL ChkMemErr('MAIN','T',IZERO)
ALLOCATE(DT_SYNC(NMESHES),STAT=IZERO)
CALL ChkMemErr('MAIN','DT_SYNC',IZERO)
ALLOCATE(DT_NEXT_SYNC(NMESHES),STAT=IZERO)
CALL ChkMemErr('MAIN','DT_NEXT_SYNC',IZERO)
ALLOCATE(MESH_STOP_STATUS(NMESHES),STAT=IZERO)
CALL ChkMemErr('MAIN','MESH_STOP_STATUS',IZERO)
! Set up dummy arrays to hold various arrays that must be exchanged among meshes
ALLOCATE(TI_LOC(N_DEVC),STAT=IZERO)
CALL ChkMemErr('MAIN','TI_LOC',IZERO)
ALLOCATE(TI_GLB(N_DEVC),STAT=IZERO)
CALL ChkMemErr('MAIN','TI_GLB',IZERO)
ALLOCATE(STATE_GLB(N_DEVC),STAT=IZERO)
CALL ChkMemErr('MAIN','STATE_GLB',IZERO)
ALLOCATE(STATE_LOC(N_DEVC),STAT=IZERO)
CALL ChkMemErr('MAIN','STATE_LOC',IZERO)
ALLOCATE(TC_GLB(N_DEVC),STAT=IZERO)
CALL ChkMemErr('MAIN','TC_GLB',IZERO)
ALLOCATE(TC_LOC(N_DEVC),STAT=IZERO)
CALL ChkMemErr('MAIN','TC_LOC',IZERO)
! Allocate a few arrays needed to exchange divergence and pressure info among meshes
IF (N_ZONE > 0) THEN
ALLOCATE(DSUM_ALL(N_ZONE),STAT=IZERO)
ALLOCATE(PSUM_ALL(N_ZONE),STAT=IZERO)
ALLOCATE(USUM_ALL(N_ZONE),STAT=IZERO)
ALLOCATE(CONNECTED_ZONES_GLOBAL(0:N_ZONE,0:N_ZONE),STAT=IZERO)
ALLOCATE(DSUM_ALL_LOCAL(N_ZONE),STAT=IZERO)
ALLOCATE(PSUM_ALL_LOCAL(N_ZONE),STAT=IZERO)
ALLOCATE(USUM_ALL_LOCAL(N_ZONE),STAT=IZERO)
ALLOCATE(CONNECTED_ZONES_LOCAL(0:N_ZONE,0:N_ZONE),STAT=IZERO)
ENDIF
! Start the clock
T = T_BEGIN
MESH_STOP_STATUS = NO_STOP
! Use the same TAG for all non-EVAC communications
TAG = 99
! Initialize global parameters
CALL INITIALIZE_GLOBAL_VARIABLES
CALL MPI_BARRIER(MPI_COMM_WORLD, IERR)
! Initialize radiation
IF (RADIATION) CALL INIT_RADIATION
CALL MPI_BARRIER(MPI_COMM_WORLD, IERR)
! Allocate and initialize mesh-specific variables
DO NM=1,NMESHES
IF (PROCESS(NM)==MYID) CALL INITIALIZE_MESH_VARIABLES(NM)
ENDDO
IF (USE_MPI) THEN
CALL MPI_ALLREDUCE(PROCESS_STOP_STATUS,STOP_STATUS,1,MPI_INTEGER,MPI_MAX,MPI_COMM_WORLD,IERR)
ELSE
STOP_STATUS = PROCESS_STOP_STATUS
ENDIF
IF (STOP_STATUS/=NO_STOP) CALL END_FDS
CALL MPI_BARRIER(MPI_COMM_WORLD, IERR)
! Allocate and initialize OMESH arrays to hold "other mesh" data for a given mesh
N_COMMUNICATIONS = 0
DO NM=1,NMESHES
IF (PROCESS(NM)==MYID) CALL INITIALIZE_MESH_EXCHANGE(NM)
ENDDO
CALL MPI_BARRIER(MPI_COMM_WORLD, IERR)
! Allocate "request" arrays to keep track of MPI communications
ALLOCATE(REQ(N_COMMUNICATIONS*40))
ALLOCATE(PREQ(N_COMMUNICATIONS*4))
REQ = MPI_REQUEST_NULL
PREQ = MPI_REQUEST_NULL
! Exchange information related to size of OMESH arrays
DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID) CYCLE
IF (EVACUATION_ONLY(NM)) CYCLE
DO NOM=1,NMESHES
IF (EVACUATION_ONLY(NOM)) CYCLE
IF (USE_MPI .AND. NM/=NOM) THEN
CALL MPI_SEND(MESHES(NM)%OMESH(NOM)%I_MIN_R,1,MPI_INTEGER,PROCESS(NOM),1,MPI_COMM_WORLD,IERR)
CALL MPI_SEND(MESHES(NM)%OMESH(NOM)%I_MAX_R,1,MPI_INTEGER,PROCESS(NOM),1,MPI_COMM_WORLD,IERR)
CALL MPI_SEND(MESHES(NM)%OMESH(NOM)%J_MIN_R,1,MPI_INTEGER,PROCESS(NOM),1,MPI_COMM_WORLD,IERR)
CALL MPI_SEND(MESHES(NM)%OMESH(NOM)%J_MAX_R,1,MPI_INTEGER,PROCESS(NOM),1,MPI_COMM_WORLD,IERR)
CALL MPI_SEND(MESHES(NM)%OMESH(NOM)%K_MIN_R,1,MPI_INTEGER,PROCESS(NOM),1,MPI_COMM_WORLD,IERR)
CALL MPI_SEND(MESHES(NM)%OMESH(NOM)%K_MAX_R,1,MPI_INTEGER,PROCESS(NOM),1,MPI_COMM_WORLD,IERR)
CALL MPI_SEND(MESHES(NM)%OMESH(NOM)%NIC_S, 1,MPI_INTEGER,PROCESS(NOM),1,MPI_COMM_WORLD,IERR)
ELSE
MESHES(NOM)%OMESH(NM)%I_MIN_S = MESHES(NM)%OMESH(NOM)%I_MIN_R
MESHES(NOM)%OMESH(NM)%I_MAX_S = MESHES(NM)%OMESH(NOM)%I_MAX_R
MESHES(NOM)%OMESH(NM)%J_MIN_S = MESHES(NM)%OMESH(NOM)%J_MIN_R
MESHES(NOM)%OMESH(NM)%J_MAX_S = MESHES(NM)%OMESH(NOM)%J_MAX_R
MESHES(NOM)%OMESH(NM)%K_MIN_S = MESHES(NM)%OMESH(NOM)%K_MIN_R
MESHES(NOM)%OMESH(NM)%K_MAX_S = MESHES(NM)%OMESH(NOM)%K_MAX_R
MESHES(NOM)%OMESH(NM)%NIC_R = MESHES(NM)%OMESH(NOM)%NIC_S
ENDIF
ENDDO
ENDDO
DO NM=1,NMESHES
IF (EVACUATION_ONLY(NM)) CYCLE
DO NOM=1,NMESHES
IF (PROCESS(NOM)/=MYID) CYCLE
IF (EVACUATION_ONLY(NOM)) CYCLE
IF (USE_MPI .AND. NM/=NOM) THEN
CALL MPI_RECV(MESHES(NOM)%OMESH(NM)%I_MIN_S,1,MPI_INTEGER,PROCESS(NM),1,MPI_COMM_WORLD,STATUS,IERR)
CALL MPI_RECV(MESHES(NOM)%OMESH(NM)%I_MAX_S,1,MPI_INTEGER,PROCESS(NM),1,MPI_COMM_WORLD,STATUS,IERR)
CALL MPI_RECV(MESHES(NOM)%OMESH(NM)%J_MIN_S,1,MPI_INTEGER,PROCESS(NM),1,MPI_COMM_WORLD,STATUS,IERR)
CALL MPI_RECV(MESHES(NOM)%OMESH(NM)%J_MAX_S,1,MPI_INTEGER,PROCESS(NM),1,MPI_COMM_WORLD,STATUS,IERR)
CALL MPI_RECV(MESHES(NOM)%OMESH(NM)%K_MIN_S,1,MPI_INTEGER,PROCESS(NM),1,MPI_COMM_WORLD,STATUS,IERR)
CALL MPI_RECV(MESHES(NOM)%OMESH(NM)%K_MAX_S,1,MPI_INTEGER,PROCESS(NM),1,MPI_COMM_WORLD,STATUS,IERR)
CALL MPI_RECV(MESHES(NOM)%OMESH(NM)%NIC_R, 1,MPI_INTEGER,PROCESS(NM),1,MPI_COMM_WORLD,STATUS,IERR)
ENDIF
ENDDO
ENDDO
! Process unusual cases where one mesh sees another but not vis verse.
DO NM=1,NMESHES
IF (PROCESS(NM)==MYID) CALL DOUBLE_CHECK(NM)
ENDDO
CALL MPI_BARRIER(MPI_COMM_WORLD, IERR)
! Initialize Mesh Exchange Arrays (All Nodes)
N_REQ=0 ! Counter for MPI requests
CALL POST_RECEIVES(0)
CALL MESH_EXCHANGE(0)
CALL MPI_BARRIER(MPI_COMM_WORLD, IERR)
! Initialize ScaRC solver
IF (PRES_METHOD == 'SCARC') CALL SCARC_SETUP
! Initialize turb arrays
DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID) CYCLE
IF (.NOT.EVACUATION_ONLY(NM)) CALL INIT_TURB_ARRAYS(NM)
ENDDO
! Initialize unstructured geometry
IF (IMMERSED_BOUNDARY_METHOD>=0) THEN
IF (CUTCELL_TEST == 0) THEN
CALL GET_CUTCELL_AREA()
!STOP
ENDIF
DO NM=1,NMESHES
IF (EVACUATION_ONLY(NM)) CYCLE
IF (PROCESS(NM)/=MYID) CYCLE
!CALL INIT_IBM(0._EB,NM)
ENDDO
ENDIF
! Initialize the flow field with random noise to eliminate false symmetries
IF (NOISE .OR. PERIODIC_TEST>0) THEN
DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID) CYCLE
IF (NOISE) CALL INITIAL_NOISE(NM)
IF (PERIODIC_TEST==1) CALL NS_ANALYTICAL_SOLUTION(NM)
IF (PERIODIC_TEST==2) CALL UVW_INIT(NM,UVW_FILE)
IF (PERIODIC_TEST==3) CALL COMPRESSION_WAVE(NM,0._EB,3)
IF (PERIODIC_TEST==4) CALL COMPRESSION_WAVE(NM,0._EB,4)
IF (PERIODIC_TEST==5) CALL STRATIFIED_MIXING_LAYER(NM)
IF (UVW_RESTART) CALL UVW_INIT(NM,CSVFINFO(NM)%UVWFILE)
ENDDO
CALL POST_RECEIVES(6)
CALL MESH_EXCHANGE(6)
PREDICTOR = .FALSE.
CORRECTOR = .TRUE.
DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID) CYCLE
CALL MATCH_VELOCITY(NM)
CALL SYNTHETIC_EDDY_SETUP(NM)
CALL VELOCITY_BC(T_BEGIN,NM)
ENDDO
ENDIF
! Potentially read data from a previous calculation
DO NM=1,NMESHES
IF (RESTART .AND. PROCESS(NM)==MYID) CALL READ_RESTART(T(NM),NM)
ENDDO
CALL MPI_BARRIER(MPI_COMM_WORLD, IERR)
! Initialize output files containing global data (Master Node Only)
IF (MYID==0) CALL INITIALIZE_GLOBAL_DUMPS
CALL INIT_EVAC_DUMPS
CALL MPI_BARRIER(MPI_COMM_WORLD, IERR)
! Initialize output files that are mesh-specific
DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID) CYCLE
CALL INITIALIZE_MESH_DUMPS(NM)
CALL INITIALIZE_PARTICLES(NM)
CALL INSERT_PARTICLES(T_BEGIN,NM)
CALL INITIALIZE_RAISED_VEG(NM)
ENDDO
CALL DEALLOCATE_VEG_ARRAYS
CALL MPI_BARRIER(MPI_COMM_WORLD, IERR)
! Initialize EVACuation routines
IF (ANY(EVACUATION_GRID)) THEN
CALL INITIALIZE_EVAC
IF (.NOT.USE_MPI .OR. (USE_MPI .AND. MYID==EVAC_PROCESS)) CALL INIT_EVAC_GROUPS(MESH_STOP_STATUS)
ENDIF
! Write out character strings to .smv file
CALL WRITE_STRINGS
! Make an initial dump of ambient values
IF (.NOT.RESTART) THEN
DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID) CYCLE
CALL UPDATE_GLOBAL_OUTPUTS(T(NM),NM)
CALL DUMP_MESH_OUTPUTS(T(NM),NM)
ENDDO
ENDIF
! Ensure the time is known to all meshes
IF (USE_MPI) THEN
REAL_BUFFER_1 = T
CALL MPI_ALLGATHERV(REAL_BUFFER_1(DISPLS(MYID)+1),COUNTS(MYID),MPI_DOUBLE_PRECISION, &
T,COUNTS,DISPLS,MPI_DOUBLE_PRECISION,MPI_COMM_WORLD,IERR)
CALL MPI_BARRIER(MPI_COMM_WORLD, IERR)
ENDIF
! If there are zones and HVAC pass PSUM
IF (HVAC_SOLVE .AND. N_ZONE>0) CALL EXCHANGE_DIVERGENCE_INFO
! Make an initial dump of global output quantities
IF (.NOT.RESTART) THEN
CALL UPDATE_CONTROLS(T)
CALL DUMP_GLOBAL_OUTPUTS(T(1))
ENDIF
! Check for changes in VENT or OBSTruction control and device status at t=T_BEGIN
IF (.NOT.RESTART) THEN
DO NM=1,NMESHES
IF (PROCESS(NM)==MYID) CALL OPEN_AND_CLOSE(T(NM),NM)
ENDDO
ENDIF
IF (ANY(EVACUATION_GRID)) THEN
IF (USE_MPI) THEN
CALL MPI_ALLREDUCE(PROCESS_STOP_STATUS,STOP_STATUS,1,MPI_INTEGER,MPI_MAX,MPI_COMM_WORLD,IERR)
ELSE
STOP_STATUS = PROCESS_STOP_STATUS
ENDIF
IF (STOP_STATUS/=NO_STOP) CALL END_FDS
CALL MPI_BARRIER(MPI_COMM_WORLD, IERR)
ENDIF
! Start the clock for time stepping
WALL_CLOCK_START_ITERATIONS = WALL_CLOCK_TIME()
!***********************************************************************************************************************************
! MAIN TIMESTEPPING LOOP
!***********************************************************************************************************************************
! General initialization of Level Set for firespread and fire front propagation without use of CFD.
IF (ACTIVE_MESH(1) .AND. PROCESS(1)==MYID) THEN
IF (VEG_LEVEL_SET_UNCOUPLED .OR. VEG_LEVEL_SET_COUPLED) CALL INITIALIZE_LEVEL_SET_FIREFRONT(1)
IF (VEG_LEVEL_SET_UNCOUPLED) THEN
CALL LEVEL_SET_FIREFRONT_PROPAGATION(T(1),1)
CALL END_LEVEL_SET
CALL END_FDS
ENDIF
ENDIF
MAIN_LOOP: DO
ICYC = ICYC + 1 ! Time step iterations
N_REQ = 0 ! Counter for MPI requests
IF (MOD(ICYC,3)==0 .AND. TIMING) THEN
MPI_WALL_TIME = MPI_WTIME() - MPI_WALL_TIME_START
WRITE(LU_ERR,'(A,I3,A,I6,A,F12.4)') ' Process ',MYID,' starts iteration',ICYC,' at ', MPI_WALL_TIME
ENDIF
! Determine if radiation or evacuation info is to be exchanged this time step
EXCHANGE_RADIATION = .FALSE.
IF (RADIATION) THEN
IF (MOD(ICYC,TIME_STEP_INCREMENT)==0) EXCHANGE_RADIATION = .TRUE.
ENDIF
EXCHANGE_EVACUATION = .FALSE.
! Synchronize clocks
IF (USE_MPI) THEN
REAL_BUFFER_1 = T
CALL MPI_ALLGATHERV(REAL_BUFFER_1(DISPLS(MYID)+1),COUNTS(MYID),MPI_DOUBLE_PRECISION, &
T,COUNTS,DISPLS,MPI_DOUBLE_PRECISION,MPI_COMM_WORLD,IERR)
ENDIF
! Check for program stops
CALL MPI_BARRIER(MPI_COMM_WORLD, IERR)
INQUIRE(FILE=FN_STOP,EXIST=EX)
IF (USE_MPI) THEN
CALL MPI_ALLREDUCE(EX,EX_GLOBAL,1,MPI_LOGICAL,MPI_LOR,MPI_COMM_WORLD,IERR)
ELSE
EX_GLOBAL = EX
ENDIF
EX = EX_GLOBAL
IF (EX) MESH_STOP_STATUS = USER_STOP
IF (EX .AND. USE_MPI) WRITE(LU_ERR,'(A,A)') PNAME(1:PNAMELEN),' read the stop file'
IF (EX .AND. .NOT.USE_MPI) WRITE(LU_ERR,'(A)') ' Stop file detected'
! Figure out fastest and slowest meshes
T_MAX = MAXVAL(T,MASK=.NOT.EVACUATION_ONLY)
T_MIN = MINVAL(T,MASK=.NOT.EVACUATION_ONLY)
IF (ALL(EVACUATION_ONLY)) T_MAX = T_EVAC
IF (ALL(EVACUATION_ONLY)) T_MIN = T_EVAC
DO NM=1,NMESHES
IF (MESH_STOP_STATUS(NM)/=NO_STOP) PROCESS_STOP_STATUS = MESH_STOP_STATUS(NM)
ENDDO
! Determine time step
IF (SYNCHRONIZE) THEN
DO NM=1,NMESHES
IF (PROCESS(NM)==MYID) DT_SYNC(NM) = MESHES(NM)%DT_NEXT
ENDDO
IF (USE_MPI) THEN
REAL_BUFFER_1 = DT_SYNC
CALL MPI_ALLGATHERV(REAL_BUFFER_1(DISPLS(MYID)+1),COUNTS(MYID),MPI_DOUBLE_PRECISION,DT_SYNC,COUNTS,DISPLS, &
MPI_DOUBLE_PRECISION,MPI_COMM_WORLD,IERR)
ENDIF
DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID) CYCLE
IF (SYNC_TIME_STEP(NM)) THEN
MESHES(NM)%DT_NEXT = MINVAL(DT_SYNC,MASK=SYNC_TIME_STEP)
T(NM) = MINVAL(T,MASK=SYNC_TIME_STEP)
ACTIVE_MESH(NM) = .TRUE.
ELSE
ACTIVE_MESH(NM) = .FALSE.
IF (T(NM)+MESHES(NM)%DT_NEXT <= T_MAX) ACTIVE_MESH(NM) = .TRUE.
IF (PROCESS_STOP_STATUS/=NO_STOP) ACTIVE_MESH(NM) = .TRUE.
ENDIF
ENDDO
ELSE
ACTIVE_MESH = .FALSE.
DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID) CYCLE
IF (T(NM)+MESHES(NM)%DT_NEXT <= T_MAX) ACTIVE_MESH(NM) = .TRUE.
IF (PROCESS_STOP_STATUS/=NO_STOP) ACTIVE_MESH(NM) = .TRUE.
ENDDO
ENDIF
! Determine when to dump out diagnostics to the .out file
DIAGNOSTICS = .FALSE.
LO10 = LOG10(REAL(MAX(1,ABS(ICYC)),EB))
IF (MOD(ICYC,10**LO10)==0 .OR. MOD(ICYC,100)==0 .OR. T_MIN>=T_END .OR. PROCESS_STOP_STATUS/=NO_STOP) DIAGNOSTICS = .TRUE.
! Give every processor the full ACTIVE_MESH array
IF (USE_MPI) THEN
LOGICAL_BUFFER_1 = ACTIVE_MESH
CALL MPI_ALLGATHERV(LOGICAL_BUFFER_1(DISPLS(MYID)+1),COUNTS(MYID),MPI_LOGICAL, &
ACTIVE_MESH,COUNTS,DISPLS,MPI_LOGICAL,MPI_COMM_WORLD,IERR)
ENDIF
! If no meshes are due to be updated, update them all
IF (ALL(.NOT.ACTIVE_MESH)) ACTIVE_MESH = .TRUE.
! If evacuation, set up special time iteration parameters
CALL EVAC_MAIN_LOOP
!============================================================================================================================
! Start of Predictor part of time step
!============================================================================================================================
PREDICTOR = .TRUE.
CORRECTOR = .FALSE.
! Diagnostic calls
DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID) CYCLE
IF (DEBUG) WRITE(LU_ERR,'(A,I8,A,I3,A,L2)') 'Cycle ',ICYC,', Mesh ',NM, ' status=',ACTIVE_MESH(NM)
IF (MOD(ICYC,3)==0 .AND. TIMING .AND. ACTIVE_MESH(NM)) THEN
MPI_WALL_TIME = MPI_WTIME() - MPI_WALL_TIME_START
WRITE(LU_ERR,'(A,I3,A,F12.4)') ' Mesh ',NM,' is active at ', MPI_WALL_TIME
ENDIF
ENDDO
! Begin the finite differencing of the PREDICTOR step
COMPUTE_FINITE_DIFFERENCES_1: DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID .OR. .NOT.ACTIVE_MESH(NM)) CYCLE COMPUTE_FINITE_DIFFERENCES_1
MESHES(NM)%DT = MESHES(NM)%DT_NEXT
NTCYC(NM) = NTCYC(NM) + 1
CALL INSERT_PARTICLES(T(NM),NM)
CALL COMPUTE_VELOCITY_FLUX(T(NM),NM,1)
CALL MASS_FINITE_DIFFERENCES(NM)
ENDDO COMPUTE_FINITE_DIFFERENCES_1
! Estimate quantities at next time step, and decrease/increase time step if necessary based on CFL condition
FIRST_PASS = .TRUE.
CHANGE_TIME_STEP_LOOP: DO
IF (FIRST_PASS .OR. SYNCHRONIZE) CALL POST_RECEIVES(1)
! Predict density and mass fractions at next time step, and then start the divergence calculation
COMPUTE_DENSITY_LOOP: DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID .OR. .NOT.ACTIVE_MESH(NM)) CYCLE COMPUTE_DENSITY_LOOP
CALL DENSITY(NM)
ENDDO COMPUTE_DENSITY_LOOP
! Exchange density and species mass fractions in interpolated boundaries
IF (FIRST_PASS .OR. SYNCHRONIZE) CALL MESH_EXCHANGE(1)
! Do mass and energy boundary conditions, and begin divergence calculation
COMPUTE_DIVERGENCE_LOOP: DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID .OR. .NOT.ACTIVE_MESH(NM)) CYCLE COMPUTE_DIVERGENCE_LOOP
IF (IMMERSED_BOUNDARY_METHOD>=0 .AND. FIRST_PASS) CALL INIT_IBM(T(NM),NM)
CALL COMPUTE_VELOCITY_FLUX(T(NM),NM,2)
CALL UPDATE_PARTICLES(T(NM),NM)
IF (FIRST_PASS .AND. HVAC_SOLVE) CALL HVAC_BC_IN(NM)
ENDDO COMPUTE_DIVERGENCE_LOOP
IF (HVAC_SOLVE) THEN
IF (FIRST_PASS .AND. USE_MPI) CALL EXCHANGE_HVAC_BC
IF (PROCESS(1)==MYID ) CALL HVAC_CALC(T(1))
IF (USE_MPI) CALL EXCHANGE_HVAC_SOLUTION
ENDIF
COMPUTE_WALL_BC_LOOP_A: DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID .OR. .NOT.ACTIVE_MESH(NM)) CYCLE COMPUTE_WALL_BC_LOOP_A
CALL WALL_BC(T(NM),NM)
CALL DIVERGENCE_PART_1(T(NM),NM)
ENDDO COMPUTE_WALL_BC_LOOP_A
! If there are pressure ZONEs, exchange integrated quantities mesh to mesh for use in the divergence calculation
IF (N_ZONE>0) CALL EXCHANGE_DIVERGENCE_INFO
! Finish the divergence calculation
FINISH_DIVERGENCE_LOOP: DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID .OR. .NOT.ACTIVE_MESH(NM)) CYCLE FINISH_DIVERGENCE_LOOP
CALL DIVERGENCE_PART_2(NM)
ENDDO FINISH_DIVERGENCE_LOOP
! Solve for the pressure at the current time step
CALL PRESSURE_ITERATION_SCHEME
CALL EVAC_PRESSURE_ITERATION_SCHEME
! Predict the velocity components at the next time step
PREDICT_VELOCITY_LOOP: DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID .OR. .NOT.ACTIVE_MESH(NM)) CYCLE PREDICT_VELOCITY_LOOP
CALL VELOCITY_PREDICTOR(T(NM)+MESHES(NM)%DT,NM,MESH_STOP_STATUS(NM))
ENDDO PREDICT_VELOCITY_LOOP
! Exchange information about the time step status, and if need be, repeat the CHANGE_TIME_STEP_LOOP
SYNCHRONIZE_ONLY: IF (SYNCHRONIZE) THEN
DISP = DISPLS(MYID)+1
IF (USE_MPI) THEN
LOGICAL_BUFFER_1 = CHANGE_TIME_STEP
CALL MPI_ALLGATHERV(LOGICAL_BUFFER_1(DISP),COUNTS(MYID),MPI_LOGICAL,CHANGE_TIME_STEP,COUNTS,DISPLS, &
MPI_LOGICAL,MPI_COMM_WORLD,IERR)
INTEGER_BUFFER_1 = MESH_STOP_STATUS
CALL MPI_ALLGATHERV(INTEGER_BUFFER_1(DISP),COUNTS(MYID),MPI_INTEGER,MESH_STOP_STATUS,COUNTS,DISPLS, &
MPI_INTEGER,MPI_COMM_WORLD,IERR)
ENDIF
IF (ANY(MESH_STOP_STATUS/=NO_STOP)) THEN
IF (ANY(MESH_STOP_STATUS==INSTABILITY_STOP)) THEN
PROCESS_STOP_STATUS = INSTABILITY_STOP
DIAGNOSTICS = .TRUE.
ENDIF
EXIT CHANGE_TIME_STEP_LOOP
ENDIF
IF (ANY(CHANGE_TIME_STEP)) THEN
CHANGE_TIME_STEP = .TRUE.
DO NM=1,NMESHES
IF (EVACUATION_ONLY(NM)) CHANGE_TIME_STEP(NM) = .FALSE.
IF (PROCESS(NM)/=MYID) CYCLE
DT_SYNC(NM) = MESHES(NM)%DT
DT_NEXT_SYNC(NM) = MESHES(NM)%DT_NEXT
ENDDO
IF (USE_MPI) THEN
REAL_BUFFER_1 = DT_SYNC
CALL MPI_ALLGATHERV(REAL_BUFFER_1(DISP),COUNTS(MYID),MPI_DOUBLE_PRECISION,DT_SYNC,COUNTS,DISPLS, &
MPI_DOUBLE_PRECISION,MPI_COMM_WORLD,IERR)
REAL_BUFFER_1 = DT_NEXT_SYNC
CALL MPI_ALLGATHERV(REAL_BUFFER_1(DISP),COUNTS(MYID),MPI_DOUBLE_PRECISION,DT_NEXT_SYNC,COUNTS,DISPLS, &
MPI_DOUBLE_PRECISION,MPI_COMM_WORLD,IERR)
ENDIF
DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID) CYCLE
IF (EVACUATION_ONLY(NM)) CYCLE
MESHES(NM)%DT_NEXT = MINVAL(DT_NEXT_SYNC,MASK=SYNC_TIME_STEP)
MESHES(NM)%DT = MINVAL(DT_SYNC,MASK=SYNC_TIME_STEP)
ENDDO
ENDIF
ENDIF SYNCHRONIZE_ONLY
DO NM=1,NMESHES
IF (PROCESS(NM)==MYID .AND. MESH_STOP_STATUS(NM)/=NO_STOP) THEN
PROCESS_STOP_STATUS = MESH_STOP_STATUS(NM)
IF (PROCESS_STOP_STATUS==INSTABILITY_STOP) DIAGNOSTICS = .TRUE.
EXIT CHANGE_TIME_STEP_LOOP
ENDIF
ENDDO
IF (.NOT.ANY(CHANGE_TIME_STEP)) EXIT CHANGE_TIME_STEP_LOOP
FIRST_PASS = .FALSE.
ENDDO CHANGE_TIME_STEP_LOOP
CHANGE_TIME_STEP = .FALSE.
! Exchange velocity and pressures at interpolated boundaries
CALL POST_RECEIVES(3)
CALL MESH_EXCHANGE(3)
! Force normal components of velocity to match at interpolated boundaries
DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID .OR. .NOT.ACTIVE_MESH(NM) .OR. MESHES(NM)%MESH_LEVEL/=0) CYCLE
CALL MATCH_VELOCITY(NM)
ENDDO
! Apply tangential velocity boundary conditions
VELOCITY_BC_LOOP: DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID .OR. .NOT.ACTIVE_MESH(NM)) CYCLE VELOCITY_BC_LOOP
CALL SYNTHETIC_TURBULENCE(MESHES(NM)%DT,T(NM),NM)
CALL VELOCITY_BC(T(NM),NM)
ENDDO VELOCITY_BC_LOOP
! Advance the time to start the CORRECTOR step
UPDATE_TIME: DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID .OR. .NOT.ACTIVE_MESH(NM)) CYCLE UPDATE_TIME
T(NM) = T(NM) + MESHES(NM)%DT
ENDDO UPDATE_TIME
! Ensure the time is known to all meshes
IF (USE_MPI) THEN
REAL_BUFFER_1 = T
CALL MPI_ALLGATHERV(REAL_BUFFER_1(DISPLS(MYID)+1),COUNTS(MYID),MPI_DOUBLE_PRECISION, &
T,COUNTS,DISPLS,MPI_DOUBLE_PRECISION,MPI_COMM_WORLD,IERR)
ENDIF
T_MAX = MAXVAL(T,MASK=.NOT.EVACUATION_ONLY)
T_MIN = MINVAL(T,MASK=.NOT.EVACUATION_ONLY)
IF (ALL(EVACUATION_ONLY)) T_MAX = T_EVAC
IF (ALL(EVACUATION_ONLY)) T_MIN = T_EVAC
!============================================================================================================================
! Start of Corrector part of time step
!============================================================================================================================
CORRECTOR = .TRUE.
PREDICTOR = .FALSE.
CALL POST_RECEIVES(4)
! Finite differences for mass and momentum equations for the second half of the time step
COMPUTE_FINITE_DIFFERENCES_2: DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID) CYCLE COMPUTE_FINITE_DIFFERENCES_2
CALL OPEN_AND_CLOSE(T(NM),NM)
IF (.NOT.ACTIVE_MESH(NM)) CYCLE COMPUTE_FINITE_DIFFERENCES_2
CALL COMPUTE_VELOCITY_FLUX(T(NM),NM,1)
CALL MASS_FINITE_DIFFERENCES(NM)
CALL DENSITY(NM)
ENDDO COMPUTE_FINITE_DIFFERENCES_2
! Exchange density and mass species
CALL MESH_EXCHANGE(4)
! Apply mass and species boundary conditions, update radiation, particles, and re-compute divergence
COMPUTE_DIVERGENCE_2: DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID .OR. .NOT.ACTIVE_MESH(NM)) CYCLE COMPUTE_DIVERGENCE_2
CALL COMPUTE_VELOCITY_FLUX(T(NM),NM,2)
CALL UPDATE_PARTICLES(T(NM),NM)
IF(.NOT. WIND_ONLY) CALL RAISED_VEG_MASS_ENERGY_TRANSFER(T(NM),NM)
! IF (HVAC_SOLVE) CALL HVAC_BC_IN(NM)
IF (N_REACTIONS > 0) CALL COMBUSTION (NM)
ENDDO COMPUTE_DIVERGENCE_2
IF (HVAC_SOLVE) THEN
! IF (USE_MPI) CALL EXCHANGE_HVAC_BC
IF (ACTIVE_MESH(1)) CALL HVAC_CALC(T(1))
! IF (USE_MPI) CALL EXCHANGE_HVAC_SOLUTION
ENDIF
COMPUTE_WALL_BC_2A: DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID .OR. .NOT.ACTIVE_MESH(NM)) CYCLE COMPUTE_WALL_BC_2A
CALL WALL_BC(T(NM),NM)
CALL BNDRY_VEG_MASS_ENERGY_TRANSFER(T(NM),NM)
CALL LEVEL_SET_FIREFRONT_PROPAGATION(T(NM),NM)
CALL COMPUTE_RADIATION(T(NM),NM)
CALL DIVERGENCE_PART_1(T(NM),NM)
ENDDO COMPUTE_WALL_BC_2A
! FREEZE_FIRE: DO NM=1,NMESHES
! IF (PROCESS(NM)/=MYID .OR. .NOT.ACTIVE_MESH(NM)) CYCLE FREEZE_FIRE
!! IF (T(NM) > 20._EB) THEN
! T(NM) = T(NM) + MESHES(NM)%DT
! CALL BNDRY_VEG_MASS_ENERGY_TRANSFER(T(NM),NM)
!! ENDIF
! ENDDO FREEZE_FIRE
! Exchange global pressure zone information
IF (N_ZONE>0) CALL EXCHANGE_DIVERGENCE_INFO
! Finish computing the divergence
FINISH_DIVERGENCE_LOOP_2: DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID .OR. .NOT.ACTIVE_MESH(NM)) CYCLE FINISH_DIVERGENCE_LOOP_2
CALL DIVERGENCE_PART_2(NM)
ENDDO FINISH_DIVERGENCE_LOOP_2
! Solve the pressure equation
CALL PRESSURE_ITERATION_SCHEME
CALL EVAC_PRESSURE_ITERATION_SCHEME
! Set up the last big exchange of info
CALL EVAC_MESH_EXCHANGE(T_EVAC,T_EVAC_SAVE,I_EVAC,ICYC,EXCHANGE_EVACUATION,0)
CALL POST_RECEIVES(6)
! Correct the velocity
CORRECT_VELOCITY_LOOP: DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID) CYCLE CORRECT_VELOCITY_LOOP
IF (.NOT.ACTIVE_MESH(NM)) CYCLE CORRECT_VELOCITY_LOOP
CALL VELOCITY_CORRECTOR(T(NM),NM)
IF (DIAGNOSTICS) CALL CHECK_DIVERGENCE(NM)
ENDDO CORRECT_VELOCITY_LOOP
! Exchange velocity and pressure at interpolated boundaries
CALL MESH_EXCHANGE(6)
! Force normal components of velocity to match at interpolated boundaries
DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID .OR. .NOT.ACTIVE_MESH(NM) .OR. MESHES(NM)%MESH_LEVEL/=0) CYCLE
IF (EVACUATION_ONLY(NM)) CYCLE
CALL MATCH_VELOCITY(NM)
ENDDO
! Apply velocity boundary conditions, and update values of HRR, DEVC, etc.
VELOCITY_BC_LOOP_2: DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID .OR. .NOT.ACTIVE_MESH(NM)) CYCLE VELOCITY_BC_LOOP_2
CALL VELOCITY_BC(T(NM),NM)
CALL UPDATE_GLOBAL_OUTPUTS(T(NM),NM)
ENDDO VELOCITY_BC_LOOP_2
! Check for dumping end of timestep outputs
CALL UPDATE_CONTROLS(T)
DO NM=1,NMESHES
IF (PROCESS(NM)==MYID .AND. ACTIVE_MESH(NM)) CALL DUMP_MESH_OUTPUTS(T(NM),NM)
ENDDO
CALL DUMP_GLOBAL_OUTPUTS(T_MIN)
! Exchange EVAC information among meshes
CALL EVAC_EXCHANGE
! Dump out diagnostics
IF (DIAGNOSTICS .OR. T_MIN>=T_END) THEN
CALL WRITE_STRINGS
CALL EXCHANGE_DIAGNOSTICS
IF (MYID==0) CALL WRITE_DIAGNOSTICS(T)
ENDIF
! Flush output file buffers
IF (T_MIN>=FLUSH_CLOCK .AND. FLUSH_FILE_BUFFERS) THEN
IF (MYID==0) CALL FLUSH_GLOBAL_BUFFERS
IF (MYID==MAX(0,EVAC_PROCESS)) CALL FLUSH_EVACUATION_BUFFERS
DO NM=1,NMESHES
IF (PROCESS(NM)==MYID) CALL FLUSH_LOCAL_BUFFERS(NM)
ENDDO
FLUSH_CLOCK = FLUSH_CLOCK + DT_FLUSH
ENDIF
! Stop the run
IF (T_MIN>=T_END .OR. PROCESS_STOP_STATUS/=NO_STOP) EXIT MAIN_LOOP
! Diagnostic print out
DO NM=1,NMESHES
IF (PROCESS(NM)/=MYID) CYCLE
IF (MOD(ICYC,3) ==0 .AND. TIMING) THEN
MPI_WALL_TIME = MPI_WTIME() - MPI_WALL_TIME_START
WRITE(LU_ERR,'(A,I3,A,I6,A,F12.4)') ' Mesh ',NM,' ends Iteration ',ICYC,' at ', MPI_WALL_TIME
ENDIF
ENDDO
ENDDO MAIN_LOOP
!***********************************************************************************************************************************
! END OF TIME STEPPING LOOP
!***********************************************************************************************************************************
! Gather up timings for the meshes
DO NM=1,NMESHES
IF (PROCESS(NM)==MYID) TUSED(1,NM) = SECOND() - TUSED(1,NM)
ENDDO
IF (USE_MPI) THEN
REAL_BUFFER_3 = TUSED
CALL MPI_GATHERV(REAL_BUFFER_3(1,DISPLS(MYID)+1),COUNTS_TIMERS(MYID),MPI_DOUBLE_PRECISION, &
TUSED,COUNTS_TIMERS,DISPLS_TIMERS,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,IERR)
ENDIF
IF (MYID==0) CALL TIMINGS
IF (PRES_METHOD == 'SCARC') CALL SCARC_TIMINGS
! Stop the calculation
CALL END_LEVEL_SET
CALL END_FDS
! The list of subroutines called from the main program follows
CONTAINS
SUBROUTINE PRESSURE_ITERATION_SCHEME
! Iterate calls to pressure solver until velocity tolerance is satisfied
IF (PREDICTOR) PRESSURE_ITERATIONS = 0
IF (PREDICTOR .AND. ITERATE_PRESSURE) THEN
CALL POST_RECEIVES(6)
CALL MESH_EXCHANGE(6)
ENDIF
IF (CORRECTOR .AND. ITERATE_PRESSURE) THEN
CALL POST_RECEIVES(3)