Context Navigation

module_generic.f90 @ 2424

Last change on this file since 2424 was 2359, checked in by lfita, 6 years ago
Adding `Flike' in different subroutines to take into account indices from non-Fortran origins
File size: 34.0 KB

Rev	Line
[715]	1	MODULE module_generic
	2	! Module with generic functions
[1608]	3
	4	!!!!!!! Subroutines/Functions
[2340]	5	! continguos_homogene_zones: Subroutine to look for contiguous zones by looking by continuous grid points
[1655]	6	! freeunit: provides the number of a free unit in which open a file
[2357]	7	! from_coordlist_2DRKmatrix: Subroutine to construct a 2D RK matrix from a list of values accompaigned
	8	! by a list of coordinates to find i,j grid-point coordinates by minimum distance
[2356]	9	! from_ptlist_2DRKmatrix: Subroutine to construct a 2D RK matrix from a list of values accompaigned
[2357]	10	! by a list of grid-point coordinates
	11	! from_ptlist_2DRKNmatrix: Subroutine to construct N 2D RK matrix from a list of values accompaigned
	12	! by a list of grid-point coordinates
[1655]	13	! GetInNamelist: Subroutine to get a paramter from a namelistfile
[2340]	14	! get_xyconlimits: Subroutine for getting the limits of contiguous values from a given point in a 2D matrix
[1612]	15	! index_list_coordsI: Function to provide the index of a given coordinate within a list of integer coordinates
[1615]	16	! Index1DArrayI: Function to provide the first index of a given value inside a 1D integer array
[1313]	17	! Index1DArrayR: Function to provide the first index of a given value inside a 1D real array
[1608]	18	! Index1DArrayR_K: Function to provide the first index of a given value inside a 1D real(r_k) array
[2340]	19	! Index1DArrayL: Function to provide the first index of a given value inside a 1D boolean array
[1313]	20	! Index2DArrayR: Function to provide the first index of a given value inside a 2D real array
[1653]	21	! Index2DArrayR_K: Function to provide the first index of a given value inside a 2D real(r_k) array
[1660]	22	! Nvalues_2DArrayI: Number of different values of a 2D integer array
[1619]	23	! mat2DPosition: Function to provide the i, j indices of a given value inside a 2D matrix
[2340]	24	! multi_Index1DArrayL: Subroutine to provide the indices of a given value inside a 1D boolean array
[2235]	25	! numberTimes: Function to provide the number of times that a given set of characters happen within a string
[1608]	26	! RangeI: Function to provide a range of d1 values from 'iniv' to 'endv', of integer values in a vector
	27	! RangeR: Function to provide a range of d1 values from 'iniv' to 'endv', of real values in a vector
	28	! RangeR_K: Function to provide a range of d1 from 'iniv' to 'endv', of real(r_k) values in a vector
[1655]	29	! stoprun: Subroutine to stop running and print a message
[2331]	30	! vectorI_S: Function to transform a vector of integers to a string of characters
	31	! vectorR_S: Function to transform a vector of reals to a string of characters
[2330]	32	! xzones_homogenization: Subroutine to homogenize 2D contiguous zones along x-axis. Zones might be
	33	! contiguous, but with different number assigned !
	34	! yzones_homogenization: Subroutine to homogenize 2D contiguous zones along y-axis. Zones might be
	35	! contiguous, but with different number assigned !
[715]	36
[1608]	37	USE module_definitions
[1655]	38	USE module_basic
[1184]	39
[715]	40	CONTAINS
	41
[1660]	42	SUBROUTINE Nvalues_2DArrayI(dx, dy, dxy, mat2DI, Nvals, vals)
	43	! Subroutine to give the number of different values of a 2D integer array
	44
	45	IMPLICIT NONE
	46
	47	INTEGER, INTENT(in) :: dx, dy, dxy
	48	INTEGER, DIMENSION(dx,dy), INTENT(in) :: mat2DI
	49	INTEGER, INTENT(out) :: Nvals
	50	INTEGER, DIMENSION(dxy), INTENT(out) :: vals
	51
	52	! Local
	53	INTEGER :: i, j, ij
	54
	55	!!!!!!! Variables
	56	! dx, dy: size of the 2D space
	57	! mat2DI: 2D integer matrix
	58	! Nvals: number of different values
	59	! vals: vector with the different values
	60
	61	fname = 'Nvalues_2DArrayI'
	62
	63	vals = 0
	64
	65	Nvals = 1
	66	vals(1) = mat2DI(1,1)
	67	DO i=1,dx
	68	DO j=1,dy
	69	IF (Index1DArrayI(vals, Nvals, mat2DI(i,j)) == -1) THEN
	70	Nvals = Nvals + 1
	71	vals(Nvals) = mat2DI(i,j)
	72	END IF
	73	END DO
	74	END DO
	75
	76	RETURN
	77
	78	END SUBROUTINE Nvalues_2DArrayI
	79
[1612]	80	INTEGER FUNCTION index_list_coordsI(Ncoords, coords, icoord)
	81	! Function to provide the index of a given coordinate within a list of integer coordinates
	82
	83	IMPLICIT NONE
	84
	85	INTEGER, INTENT(in) :: Ncoords
	86	INTEGER, DIMENSION(Ncoords,2), INTENT(in) :: coords
	87	INTEGER, DIMENSION(2), INTENT(in) :: icoord
	88
	89	! Local
	90	INTEGER, DIMENSION(Ncoords) :: dist
	91	INTEGER :: i,mindist
	92	INTEGER, DIMENSION(1) :: iloc
	93
	94	!!!!!!! Variables
	95	! Ncoords: number of coordinates in the list
	96	! coords: list of coordinates
	97	! icoord: coordinate to find
	98
	99	fname = 'index_list_coordsI'
	100
	101	dist = (coords(:,1)-icoord(1))2+(coords(:,2)-icoord(2))2
	102
	103	IF (ANY(dist == 0)) THEN
	104	iloc = MINLOC(dist)
	105	index_list_coordsI = iloc(1)
	106	ELSE
	107	index_list_coordsI = -1
	108	END IF
	109
	110	END FUNCTION index_list_coordsI
	111
[1615]	112	INTEGER FUNCTION Index1DArrayI(array1D, d1, val)
	113	! Function to provide the first index of a given value inside a 1D integer array
	114
	115	IMPLICIT NONE
	116
	117	INTEGER, INTENT(in) :: d1
	118	INTEGER, INTENT(in) :: val
	119	INTEGER, DIMENSION(d1), INTENT(in) :: array1D
	120
	121	! Local
	122	INTEGER :: i
	123
	124	fname = 'Index1DArrayI'
	125
	126	Index1DArrayI = -1
	127
	128	DO i=1,d1
	129	IF (array1d(i) == val) THEN
	130	Index1DArrayI = i
	131	EXIT
	132	END IF
	133	END DO
	134
	135	END FUNCTION Index1DArrayI
	136
[764]	137	INTEGER FUNCTION Index1DArrayR(array1D, d1, val)
	138	! Function to provide the first index of a given value inside a 1D real array
	139
	140	IMPLICIT NONE
	141
	142	INTEGER, INTENT(in) :: d1
[1653]	143	REAL, INTENT(in) :: val
	144	REAL, DIMENSION(d1), INTENT(in) :: array1D
[764]	145
	146	! Local
	147	INTEGER :: i
	148
	149	fname = 'Index1DArrayR'
	150
	151	Index1DArrayR = -1
	152
	153	DO i=1,d1
	154	IF (array1d(i) == val) THEN
	155	Index1DArrayR = i
	156	EXIT
	157	END IF
	158	END DO
	159
	160	END FUNCTION Index1DArrayR
	161
[1608]	162	INTEGER FUNCTION Index1DArrayR_K(array1D, d1, val)
	163	! Function to provide the first index of a given value inside a 1D real(r_k) array
	164
	165	IMPLICIT NONE
	166
	167	INTEGER, INTENT(in) :: d1
	168	REAL(r_k), INTENT(in) :: val
	169	REAL(r_k), DIMENSION(d1), INTENT(in) :: array1D
	170
	171	! Local
	172	INTEGER :: i
	173
	174	fname = 'Index1DArrayR_K'
	175
	176	Index1DArrayR_K = -1
	177
	178	DO i=1,d1
	179	IF (array1d(i) == val) THEN
	180	Index1DArrayR_K = i
	181	EXIT
	182	END IF
	183	END DO
	184
	185	END FUNCTION Index1DArrayR_K
	186
[2340]	187	INTEGER FUNCTION Index1DArrayL(array1D, d1, val)
	188	! Function to provide the first index of a given value inside a 1D boolean array
	189
	190	IMPLICIT NONE
	191
	192	INTEGER, INTENT(in) :: d1
	193	LOGICAL, INTENT(in) :: val
	194	LOGICAL, DIMENSION(d1), INTENT(in) :: array1D
	195
	196	! Local
	197	INTEGER :: i
	198	CHARACTER(LEN=50) :: fname
	199
	200	fname = 'Index1DArrayL'
	201
	202	Index1DArrayL = -1
	203
	204	DO i=1,d1
	205	IF (array1d(i) .EQV. val) THEN
	206	Index1DArrayL = i
	207	EXIT
	208	END IF
	209	END DO
	210
	211	END FUNCTION Index1DArrayL
	212
	213	SUBROUTINE multi_Index1DArrayL(array1D, d1, val, Ntimes, pos)
	214	! Subroutine to provide the indices of a given value inside a 1D boolean array
	215
	216	IMPLICIT NONE
	217
	218	INTEGER, INTENT(in) :: d1
	219	LOGICAL, INTENT(in) :: val
	220	LOGICAL, DIMENSION(d1), INTENT(in) :: array1D
	221	INTEGER, INTENT(out) :: Ntimes
	222	INTEGER, DIMENSION(d1), INTENT(out) :: pos
	223
	224	! Local
	225	INTEGER :: i
	226	CHARACTER(LEN=50) :: fname
	227
	228	!!!!!!! Variables
	229	! array1D: 1D Array of values
	230	! d1: length of array
	231	! val: Value to look for
	232	! Ntimes: Number of times val is found within array1D
	233	! pos: positions of the values of val
	234
	235	fname = 'multi_Index1DArrayL'
	236
	237	Ntimes = 0
	238	pos = -1
	239
	240	DO i=1,d1
	241	IF (array1d(i) .EQV. val) THEN
	242	Ntimes = Ntimes + 1
	243	pos(Ntimes) = i
	244	END IF
	245	END DO
	246
	247	END SUBROUTINE multi_Index1DArrayL
	248
[715]	249	FUNCTION Index2DArrayR(array2D, d1, d2, val)
	250	! Function to provide the first index of a given value inside a 2D real array
	251
	252	IMPLICIT NONE
	253
	254	INTEGER, INTENT(in) :: d1, d2
[1653]	255	REAL, INTENT(in) :: val
	256	REAL, DIMENSION(d1,d2), INTENT(in) :: array2D
[715]	257	INTEGER, DIMENSION(2) :: Index2DArrayR
	258
	259	! Local
	260	INTEGER :: i, j
	261
	262	fname = 'Index2DArrayR'
	263
	264	Index2DArrayR = -1
	265
	266	DO i=1,d1
	267	DO j=1,d2
	268	IF (array2d(i,j) == val) THEN
	269	Index2DArrayR(1) = i
	270	Index2DArrayR(2) = j
	271	EXIT
	272	END IF
	273	END DO
	274	END DO
	275
	276	END FUNCTION Index2DArrayR
	277
[1653]	278	FUNCTION Index2DArrayR_K(array2D, d1, d2, val)
	279	! Function to provide the first index of a given value inside a 2D real array
	280
	281	IMPLICIT NONE
	282
	283	INTEGER, INTENT(in) :: d1, d2
	284	REAL(r_k), INTENT(in) :: val
	285	REAL(r_k), DIMENSION(d1,d2), INTENT(in) :: array2D
	286	INTEGER, DIMENSION(2) :: Index2DArrayR_K
	287
	288	! Local
	289	INTEGER :: i, j
	290
	291	fname = 'Index2DArrayR_K'
	292
	293	Index2DArrayR_K = -1
	294
	295	DO i=1,d1
	296	DO j=1,d2
	297	IF (array2d(i,j) == val) THEN
	298	Index2DArrayR_K(1) = i
	299	Index2DArrayR_K(2) = j
	300	EXIT
	301	END IF
	302	END DO
	303	END DO
	304
	305	END FUNCTION Index2DArrayR_K
	306
[2235]	307	INTEGER FUNCTION numberTimes(String, charv)
	308	! Function to provide the number of times that a given set of characters happen within a string
	309
	310	IMPLICIT NONE
	311
	312	CHARACTER(LEN=*), INTENT(IN) :: String, charv
	313
	314	! Local
	315	INTEGER :: i, Lstring, Lcharv
	316
	317	numberTimes = 0
	318
	319	Lstring = LEN_TRIM(String)
	320	Lcharv = LEN_TRIM(charv)
	321
	322	DO i=1,Lstring - Lcharv
	323	IF (String(i:i+Lcharv-1) == TRIM(charv)) numberTimes = numberTimes + 1
	324	END DO
	325
	326	END FUNCTION numberTimes
	327
	328
[1608]	329	FUNCTION RangeI(d1, iniv, endv)
	330	! Function to provide a range of d1 values from 'iniv' to 'endv', of integer values in a vector
	331
	332	IMPLICIT NONE
	333
	334	INTEGER, INTENT(in) :: d1, iniv, endv
	335	INTEGER, DIMENSION(d1) :: RangeI
	336
	337	! Local
	338	INTEGER :: i, intv
	339
	340	fname = 'RangeI'
	341
	342	intv = (endv - iniv) / (d1*1 - 1)
	343
	344	RangeI(1) = iniv
	345	DO i=2,d1
	346	RangeI(i) = RangeI(i-1) + intv
	347	END DO
	348
	349	END FUNCTION RangeI
	350
	351	FUNCTION RangeR(d1, iniv, endv)
	352	! Function to provide a range of d1 from 'iniv' to 'endv', of real values in a vector
	353
	354	IMPLICIT NONE
	355
	356	INTEGER, INTENT(in) :: d1
	357	REAL, INTENT(in) :: iniv, endv
	358	REAL, DIMENSION(d1) :: RangeR
	359
	360	! Local
	361	INTEGER :: i
	362	REAL :: intv
	363
	364	fname = 'RangeR'
	365
	366	intv = (endv - iniv) / (d1*1. - 1.)
	367
	368	RangeR(1) = iniv
	369	DO i=2,d1
	370	RangeR(i) = RangeR(i-1) + intv
	371	END DO
	372
	373	END FUNCTION RangeR
	374
	375	FUNCTION RangeR_K(d1, iniv, endv)
	376	! Function to provide a range of d1 from 'iniv' to 'endv', of real(r_k) values in a vector
	377
	378	IMPLICIT NONE
	379
	380	INTEGER, INTENT(in) :: d1
	381	REAL(r_k), INTENT(in) :: iniv, endv
	382	REAL(r_k), DIMENSION(d1) :: RangeR_K
	383
	384	! Local
	385	INTEGER :: i
	386	REAL(r_k) :: intv
	387
	388	fname = 'RangeR_K'
	389
	390	intv = (endv - iniv) / (d1*oneRK-oneRK)
	391
	392	RangeR_K(1) = iniv
	393	DO i=2,d1
	394	RangeR_K(i) = RangeR_K(i-1) + intv
	395	END DO
	396
	397	END FUNCTION RangeR_K
	398
[1655]	399	INTEGER FUNCTION freeunit()
	400	! provides the number of a free unit in which open a file
[1608]	401
[1619]	402	IMPLICIT NONE
	403
[1655]	404	LOGICAL :: is_used
[1619]	405
[1655]	406	is_used = .true.
	407	DO freeunit=10,100
[1660]	408	INQUIRE(unit=freeunit, opened=is_used)
[1655]	409	IF (.not. is_used) EXIT
	410	END DO
	411
	412	RETURN
	413
	414	END FUNCTION freeunit
	415
	416	SUBROUTINE GetInNamelist(namelistfile, param, kindparam, Ival, Rval, Lval, Sval)
	417	! Subroutine to get a paramter from a namelistfile
	418
	419	IMPLICIT NONE
	420
	421	CHARACTER(LEN=*), INTENT(IN) :: namelistfile, param
	422	CHARACTER(LEN=1), INTENT(IN) :: kindparam
	423	INTEGER, OPTIONAL, INTENT(OUT) :: Ival
	424	REAL, OPTIONAL, INTENT(OUT) :: Rval
	425	LOGICAL, OPTIONAL, INTENT(OUT) :: Lval
	426	CHARACTER(LEN=200), OPTIONAL, INTENT(OUT) :: Sval
	427
	428	! Local
	429	INTEGER :: i, funit, ios
	430	INTEGER :: Lparam, posparam
	431	LOGICAL :: is_used
	432	CHARACTER(LEN=1000) :: line, message
	433	CHARACTER(LEN=200), DIMENSION(2) :: lvals
	434	CHARACTER(LEN=200) :: pval
	435
	436	!!!!!!! Variables
	437	! namelistfile: name of the namelist file
	438	! param: parameter to get
	439	! paramkind: kind of the parameter (I: Integer, L: boolean, R: Real, S: String)
	440
	441	fname = 'GetInNamelist'
	442
	443	! Reading dimensions file and defining dimensions
	444	is_used = .true.
	445	DO funit=10,100
	446	INQUIRE(unit=funit, opened=is_used)
	447	IF (.not. is_used) EXIT
	448	END DO
	449
	450	OPEN(funit, FILE=TRIM(namelistfile), STATUS='old', FORM='formatted', IOSTAT=ios)
	451	IF ( ios /= 0 ) CALL stoprun(message, fname)
	452
	453	Lparam = LEN_TRIM(param)
	454
	455	DO i=1,10000
	456	READ(funit,"(A200)",END=100)line
	457	posparam = INDEX(TRIM(line), TRIM(param))
	458	IF (posparam /= 0) EXIT
	459
	460	END DO
	461	100 CONTINUE
	462
	463	IF (posparam == 0) THEN
	464	message = "namelist '" // TRIM(namelistfile) // "' does not have parameter '" // TRIM(param) // &
	465	"' !!"
	466	CALL stoprun(message, fname)
	467	END IF
	468
	469	CLOSE(UNIT=funit)
	470
	471	CALL split(line, '=', 2, lvals)
	472	IF (kindparam /= 'S') THEN
	473	CALL RemoveNonNum(lvals(2), pval)
	474	END IF
	475
	476	! L. Fita, LMD. October 2015
	477	! Up to now, only getting scalar values
	478	kparam: SELECT CASE (kindparam)
	479	CASE ('I')
	480	Ival = StoI(pval)
	481	! PRINT *,TRIM(param),'= ', Ival
	482	CASE ('L')
	483	Lval = StoL(pval)
	484	! PRINT *,TRIM(param),'= ', Lval
	485	CASE ('R')
	486	Rval = StoR(pval)
	487	! PRINT *,TRIM(param),'= ', Rval
	488	CASE ('S')
	489	Sval = lvals(2)
	490
	491	CASE DEFAULT
	492	message = "type of parameter '" // kindparam // "' not ready !!"
	493	CALL stoprun(message, fname)
	494
	495	END SELECT kparam
	496
	497	END SUBROUTINE GetInNamelist
	498
	499	SUBROUTINE stoprun(msg, fname)
	500	! Subroutine to stop running and print a message
	501
	502	IMPLICIT NONE
	503
	504	CHARACTER(LEN=*), INTENT(IN) :: fname
	505	CHARACTER(LEN=*), INTENT(IN) :: msg
	506
	507	! local
	508	CHARACTER(LEN=50) :: errmsg, warnmsg
	509
	510	errmsg = 'ERROR -- error -- ERROR -- error'
	511
	512	PRINT *, TRIM(errmsg)
	513	PRINT *, ' ' // TRIM(fname) // ': ' // TRIM(msg)
[1683]	514	STOP
[1655]	515
	516	END SUBROUTINE stoprun
	517
[2330]	518	SUBROUTINE xzones_homogenization(dx, dy, inzones, outzones)
	519	! Subroutine to homogenize 2D contiguous zones along x-axis, zones might be contiguous, but with
	520	! different number assigned !
[2328]	521	! Here we have a 2D matrix of integers, with contiguous integer filled zones, zero outside any zone
	522	! It might be, that within the same zone, might be lines which do not share the same integer
	523	! 0 0 0 0 0 0 0 0 0 0
	524	! 0 1 1 0 0 0 1 1 0 0
	525	! 0 2 0 0 1 == > 0 1 0 0 2
	526	! 0 1 1 0 0 0 1 1 0 0
	527
	528	IMPLICIT NONE
	529
	530	INTEGER, INTENT(in) :: dx, dy
	531	INTEGER, DIMENSION(dx,dy), INTENT(in) :: inzones
	532	INTEGER, DIMENSION(dx,dy), INTENT(out) :: outzones
	533
	534	! Local
	535	INTEGER :: i,j,k
	536	INTEGER :: Nmaxzones, TOTzones
	537	LOGICAL :: assigned
[2330]	538	INTEGER, DIMENSION(dx) :: prevline
	539	INTEGER, DIMENSION(dy) :: Nxzones
	540	INTEGER, DIMENSION(:,:,:), ALLOCATABLE :: zones
	541
	542	!!!!!!! Variables
	543	! dx, dy: Shape of the 2D space
	544	! inzones: zones to homogenize
	545	! outzones: zones homogenized
	546
	547	fname = 'xzones_homogenization'
	548
	549	! Maximum possible number of zones
	550	Nmaxzones = INT((dx/2)*(dy/2))
	551
	552	! Matrix with [i,j,Nzone,izone/ezone]
	553	IF (ALLOCATED(zones)) DEALLOCATE(zones)
	554	ALLOCATE(zones(dy,Nmaxzones,3))
	555
	556	zones = 0
	557	Nxzones = 0
	558	! Getting beginning/end of y-bands
	559	DO j=1, dy
	560	k = 0
	561	i = 1
	562	IF (inzones(i,j) /= 0) THEN
	563	k = k + 1
	564	zones(j,k,1) = i
	565	zones(j,k,3) = k
	566	END IF
	567	DO i=2, dx
	568	IF ( (inzones(i,j) /= 0) .AND. (inzones(i-1,j) == 0)) THEN
	569	k = k+1
	570	zones(j,k,1) = i
	571	zones(j,k,3) = k
	572	ELSE IF ( (inzones(i-1,j) /= 0) .AND. (inzones(i,j) == 0)) THEN
	573	zones(j,k,2) = i-1
	574	zones(j,k,3) = k
	575	END IF
	576	END DO
	577	IF (k > 0) THEN
	578	IF (zones(j,k,2) == 0) zones(j,k,2) = dx
	579	END IF
	580	Nxzones(j) = k
	581	END DO
	582
	583	! Homogenizing contigous zones
	584	outzones = 0
	585	TOTzones = 0
	586	j = 1
	587	DO k = 1, Nxzones(j)
	588	TOTzones = TOTzones + 1
	589	DO i=zones(j,k,1), zones(j,k,2)
	590	outzones(i,j) = TOTzones
	591	END DO
	592	END DO
	593
	594	DO j=2, dy
	595	prevline = outzones(:,j-1)
	596	DO k = 1, Nxzones(j)
	597	assigned = .FALSE.
	598	DO i=zones(j,k,1), zones(j,k,2)
	599	IF (prevline(i) /= 0) THEN
	600	outzones(zones(j,k,1):zones(j,k,2),j) = prevline(i)
	601	assigned = .TRUE.
	602	EXIT
	603	END IF
	604	END DO
	605	IF (.NOT.assigned) THEN
	606	TOTzones = TOTzones + 1
	607	DO i=zones(j,k,1), zones(j,k,2)
	608	outzones(i,j) = TOTzones
	609	END DO
	610	END IF
	611	END DO
	612	END DO
	613
	614	IF (ALLOCATED(zones)) DEALLOCATE(zones)
	615
	616	END SUBROUTINE xzones_homogenization
	617
	618	SUBROUTINE yzones_homogenization(dx, dy, inzones, outzones)
	619	! Subroutine to homogenize 2D contiguous zones along y-axis, zones might be contiguous, but with
	620	! different number assigned !
	621	! Here we have a 2D matrix of integers, with contiguous integer filled zones, zero outside any zone
	622	! It might be, that within the same zone, might be lines which do not share the same integer
	623	! 0 0 0 0 0 0 0 0 0 0
	624	! 0 1 1 0 0 0 1 1 0 0
	625	! 0 2 0 0 1 == > 0 1 0 0 2
	626	! 0 1 1 0 0 0 1 1 0 0
	627
	628	IMPLICIT NONE
	629
	630	INTEGER, INTENT(in) :: dx, dy
	631	INTEGER, DIMENSION(dx,dy), INTENT(in) :: inzones
	632	INTEGER, DIMENSION(dx,dy), INTENT(out) :: outzones
	633
	634	! Local
	635	INTEGER :: i,j,k
	636	INTEGER :: Nmaxzones, TOTzones
	637	LOGICAL :: assigned
[2328]	638	INTEGER, DIMENSION(dy) :: prevline
	639	INTEGER, DIMENSION(dx) :: Nyzones
	640	INTEGER, DIMENSION(:,:,:), ALLOCATABLE :: zones
	641
	642	!!!!!!! Variables
	643	! dx, dy: Shape of the 2D space
	644	! inzones: zones to homogenize
	645	! outzones: zones homogenized
	646
[2330]	647	fname = 'yzones_homogenization'
[2328]	648
	649	! Maximum possible number of zones
	650	Nmaxzones = INT((dx/2)*(dy/2))
	651
	652	! Matrix with [i,j,Nzone,izone/ezone]
	653	IF (ALLOCATED(zones)) DEALLOCATE(zones)
	654	ALLOCATE(zones(dx,Nmaxzones,3))
	655
	656	zones = 0
	657	Nyzones = 0
	658	! Getting beginning/end of y-bands
	659	DO i=1, dx
	660	k = 0
	661	j = 1
	662	IF (inzones(i,j) /= 0) THEN
	663	k = k + 1
	664	zones(i,k,1) = j
	665	zones(i,k,3) = k
	666	END IF
	667	DO j=2, dy
	668	IF ( (inzones(i,j) /= 0) .AND. (inzones(i,j-1) == 0)) THEN
	669	k = k+1
	670	zones(i,k,1) = j
	671	zones(i,k,3) = k
	672	ELSE IF ( (inzones(i,j-1) /= 0) .AND. (inzones(i,j) == 0)) THEN
	673	zones(i,k,2) = j-1
	674	zones(i,k,3) = k
	675	END IF
	676	END DO
	677	IF (k > 0) THEN
	678	IF (zones(i,k,2) == 0) zones(i,k,2) = dy
	679	END IF
	680	Nyzones(i) = k
	681	END DO
	682
	683	! Homogenizing contigous zones
	684	outzones = 0
	685	TOTzones = 0
	686	i = 1
	687	DO k = 1, Nyzones(i)
	688	TOTzones = TOTzones + 1
	689	DO j=zones(i,k,1), zones(i,k,2)
	690	outzones(i,j) = TOTzones
	691	END DO
	692	END DO
	693
	694	DO i=2, dx
	695	prevline = outzones(i-1,:)
	696	DO k = 1, Nyzones(i)
	697	assigned = .FALSE.
	698	DO j=zones(i,k,1), zones(i,k,2)
	699	IF (prevline(j) /= 0) THEN
	700	outzones(i,zones(i,k,1):zones(i,k,2)) = prevline(j)
	701	assigned = .TRUE.
	702	EXIT
	703	END IF
	704	END DO
	705	IF (.NOT.assigned) THEN
	706	TOTzones = TOTzones + 1
	707	DO j=zones(i,k,1), zones(i,k,2)
	708	outzones(i,j) = TOTzones
	709	END DO
	710	END IF
	711	END DO
	712	END DO
	713
	714	IF (ALLOCATED(zones)) DEALLOCATE(zones)
	715
[2330]	716	END SUBROUTINE yzones_homogenization
[2328]	717
[2331]	718	CHARACTER(len=1000) FUNCTION vectorI_S(d1, vector)
	719	! Function to transform a vector of integers to a string of characters
	720
	721	IMPLICIT NONE
	722
	723	INTEGER, INTENT(in) :: d1
	724	INTEGER, DIMENSION(d1), INTENT(in) :: vector
	725
	726	! Local
	727	INTEGER :: iv
	728	CHARACTER(len=50) :: IS
	729
	730	!!!!!!! Variables
	731	! d1: length of the vector
	732	! vector: values to transform
	733
	734	fname = 'vectorI_S'
	735
	736	vectorI_S = ''
	737	DO iv=1, d1
	738	WRITE(IS, '(I50)')vector(iv)
	739	IF (iv == 1) THEN
	740	vectorI_S = TRIM(IS)
	741	ELSE
	742	vectorI_S = TRIM(vectorI_S) // ', ' // TRIM(IS)
	743	END IF
	744	END DO
	745
	746	END FUNCTION vectorI_S
	747
	748	CHARACTER(len=1000) FUNCTION vectorR_S(d1, vector)
	749	! Function to transform a vector of reals to a string of characters
	750
	751	IMPLICIT NONE
	752
	753	INTEGER, INTENT(in) :: d1
	754	REAL, DIMENSION(d1), INTENT(in) :: vector
	755
	756	! Local
	757	INTEGER :: iv
	758	CHARACTER(len=50) :: RS
	759
	760	!!!!!!! Variables
	761	! d1: length of the vector
	762	! vector: values to transform
	763
	764	fname = 'vectorR_S'
	765
	766	vectorR_S = ''
	767	DO iv=1, d1
	768	WRITE(RS, '(F50.25)')vector(iv)
	769	IF (iv == 1) THEN
	770	vectorR_S = TRIM(RS)
	771	ELSE
	772	vectorR_S = TRIM(vectorR_S) // ', ' // TRIM(RS)
	773	END IF
	774	END DO
	775
	776	END FUNCTION vectorR_S
	777
[2340]	778	SUBROUTINE continguos_homogene_zones(dx, dy, matvals, Nzones, contzones)
	779	! Subroutine to look for contiguous zones by looking by continuous grid points
	780
	781	IMPLICIT NONE
	782
	783	INTEGER, INTENT(in) :: dx, dy
	784	INTEGER, DIMENSION(dx,dy), INTENT(in) :: matvals
	785	INTEGER, INTENT(out) :: Nzones
	786	INTEGER, DIMENSION(dx,dy), INTENT(out) :: contzones
	787	! Local
	788	INTEGER :: i,j, k
	789	INTEGER :: ii, ei, ij, ej
	790	INTEGER :: Ncont, Nassigned, Ncontmin
	791	LOGICAL, DIMENSION(dx,dy) :: assigned, notdone
	792	INTEGER, DIMENSION(:), ALLOCATABLE :: pzones, allzones
	793	LOGICAL, DIMENSION(:), ALLOCATABLE :: passigns
	794
	795	!!!!!!! Variables
	796	! dx, dy: shape of the matrix
	797	! matvals: matrix with the values
	798	! contzones: homogeneous zones found
	799
	800	fname = 'continguos_homogene_zones'
	801
	802	! Vector to keep track of all zone values
	803	IF (ALLOCATED(allzones)) DEALLOCATE(allzones)
	804	ALLOCATE(allzones(dx*dy/4))
	805	allzones = 0
	806
	807	assigned = .FALSE.
	808	notdone = .TRUE.
	809	contzones = -1
	810	Nzones = 0
	811	DO i=1, dx
	812	DO j=1, dy
	813	! First
	814	IF (matvals(i,j) /= 0 .AND. notdone(i,j)) THEN
	815	CALL get_xyconlimits(dx, dy, matvals, 0, i, j, ii, ei, ij, ej)
	816
	817	! Has any point of the rays already been assigned?
	818	! Along x
	819	IF (ALLOCATED(pzones)) DEALLOCATE(pzones)
	820	ALLOCATE(pzones(ei-ii+1))
	821	IF (ALLOCATED(passigns)) DEALLOCATE(passigns)
	822	ALLOCATE(passigns(ei-ii+1))
	823	passigns = assigned(ii:ei,j)
	824	CALL multi_Index1DArrayL(passigns, ei-ii+1, .TRUE., Nassigned, pzones)
	825	IF (Nassigned /= 0) THEN
	826	Ncontmin = 10000000
	827	DO k=1, Nassigned
	828	IF (contzones(ii+pzones(k)-1,j) < Ncontmin) Ncontmin = contzones(ii+pzones(k)-1,j)
	829	END DO
	830	! If there is more than one assigned value change all values to the minimum one
	831	DO k=1, Nassigned
	832	IF (contzones(ii+pzones(k)-1,j) /= Ncontmin) THEN
	833	WHERE (contzones == contzones(ii+pzones(k)-1,j) .AND. contzones /= Ncontmin)
	834	contzones = Ncontmin
	835	END WHERE
	836	allzones(contzones(ii+pzones(k)-1,j)) = 0
	837	END IF
	838	END DO
	839	Ncont = Ncontmin
	840	END IF
	841	! Along y
	842	IF (ALLOCATED(pzones)) DEALLOCATE(pzones)
	843	ALLOCATE(pzones(ej-ij+1))
	844	IF (ALLOCATED(passigns)) DEALLOCATE(passigns)
	845	ALLOCATE(passigns(ej-ij+1))
	846	passigns = assigned(i,ij:ej)
	847	CALL multi_Index1DArrayL(passigns, ej-ij+1, .TRUE., Nassigned, pzones)
	848	IF (Nassigned /= 0) THEN
	849	Ncontmin = 10000000
	850	DO k=1, Nassigned
	851	IF (contzones(i,ij+pzones(k)-1) < Ncontmin) Ncontmin = contzones(i,ij+pzones(k)-1)
	852	END DO
	853	! If there is more than one assigned value change all values to the minimum one
	854	DO k=1, Nassigned
	855	IF (contzones(i,ij+pzones(k)-1) /= Ncontmin) THEN
	856	WHERE (contzones == contzones(i,ij+pzones(k)-1) .AND. contzones /= Ncontmin)
	857	contzones = Ncontmin
	858	END WHERE
	859	allzones(contzones(i,ij+pzones(k)-1)) = 0
	860	END IF
	861	END DO
	862	Ncont = Ncontmin
	863	END IF
	864
	865	IF (.NOT.assigned(i,j)) THEN
	866	Nzones = Nzones + 1
	867	Ncont = Nzones
	868	allzones(Nzones) = 1
	869	ELSE
	870	Ncont = contzones(i,j)
	871	END IF
	872	contzones(i,j) = Ncont
	873	contzones(ii:ei,j) = Ncont
	874	contzones(i,ij:ej) = Ncont
	875	notdone(i,j) = .FALSE.
	876	assigned(i,j) = .TRUE.
	877	assigned(ii:ei,j) = .TRUE.
	878	assigned(i,ij:ej) = .TRUE.
	879	END IF
	880	END DO
	881	END DO
	882
	883	! Using allzones to provide continuous assigned values
	884	Nzones = 0
	885	DO k=1, dx*dy/4
	886	IF (allzones(k) /= 0) THEN
	887	Nzones = Nzones + 1
	888	IF (allzones(k) /= Nzones) THEN
	889	WHERE(contzones == allzones(k))
	890	contzones = Nzones
	891	END WHERE
	892	END IF
	893	END IF
	894	END DO
	895
	896	RETURN
	897
	898	END SUBROUTINE continguos_homogene_zones
	899
	900	SUBROUTINE get_xyconlimits(d1, d2, matv, NOval, i, j, ix, ex, iy, ey)
	901	! Subroutine for getting the limits of contiguous values from a given point in a 2D matrix
	902
	903	IMPLICIT NONE
	904
	905	INTEGER, INTENT(in) :: d1, d2, i, j, NOval
	906	INTEGER, DIMENSION(d1,d2), INTENT(in) :: matv
	907	INTEGER, INTENT(out) :: ix, ex, iy, ey
	908
	909	! Local
	910	INTEGER :: i1, j1
	911	LOGICAL :: found
	912
	913	!!!!!!! Variables
	914	! d1, d2: Shape of input 2D values
	915	! i, j: grid point from which get the contiguous values
	916	! NOval: value for grid points without data
	917	! matv: 2D matrx values
	918	! ix, ex, iy, ey: limits from a grid point
	919
	920	fname = 'get_xyconlimits'
	921
	922	ix = -1
	923	ex = -1
	924	iy = -1
	925	ey = -1
	926
	927	! Before i
	928	IF (i > 1) THEN
	929	ix = i
	930	found = .FALSE.
	931	DO i1=i-1,1,-1
	932	IF (matv(i1,j) == NOval) THEN
	933	ix = i1 + 1
	934	found = .TRUE.
	935	EXIT
	936	END IF
	937	END DO
	938	IF (.NOT.found) ix=1
	939	ELSE
	940	ix = 1
	941	END IF
	942	! After i
	943	IF (i < d1) THEN
	944	ex = i
	945	found = .FALSE.
	946	DO i1=i+1,d1
	947	IF (matv(i1,j) == NOval) THEN
	948	ex = i1 - 1
	949	found = .TRUE.
	950	EXIT
	951	END IF
	952	END DO
	953	IF (.NOT.found) ex=d1
	954	ELSE
	955	ex = d1
	956	END IF
	957
	958	! Before j
	959	IF (j > 1) THEN
	960	iy = j
	961	found = .FALSE.
	962	DO j1=j-1,1,-1
	963	IF (matv(i,j1) == NOval) THEN
	964	iy = j1 + 1
	965	found = .TRUE.
	966	EXIT
	967	END IF
	968	END DO
	969	IF (.NOT.found) iy=1
	970	ELSE
	971	iy = 1
	972	END IF
	973	! After j
	974	IF (j < d2) THEN
	975	ey = j
	976	found = .FALSE.
	977	DO j1=j+1,d2
	978	IF (matv(i,j1) == NOval) THEN
	979	ey = j1 - 1
	980	found = .TRUE.
	981	EXIT
	982	END IF
	983	END DO
	984	IF (.NOT.found) ey=d2
	985	ELSE
	986	ey = d2
	987	END IF
	988
	989	RETURN
	990
	991	END SUBROUTINE get_xyconlimits
	992
[2359]	993	SUBROUTINE from_ptlist_2DRKmatrix(Npts, pts, Flike, vals, dx, dy, NOval, matrix)
[2357]	994	! Subroutine to construct a 2D RK matrix from a list of values accompaigned by a list of grid-point
	995	! coordinates
[2356]	996
	997	IMPLICIT NONE
	998
	999	INTEGER, INTENT(in) :: Npts, dx, dy
[2359]	1000	LOGICAL, INTENT(in) :: Flike
[2356]	1001	REAL(r_k), INTENT(in) :: NOval
	1002	INTEGER, DIMENSION(Npts,2), INTENT(in) :: pts
	1003	REAL(r_k), DIMENSION(Npts), INTENT(in) :: vals
	1004	REAL(r_k), DIMENSION(dx,dy), INTENT(out) :: matrix
	1005
	1006	! Local
	1007	INTEGER :: iv, i, j
	1008
	1009	!!!!!!! Variables
	1010	! Npts: Number of values of the list
	1011	! pts: 2D matrix coordinates of the values
[2359]	1012	! Flike: whether input is Fortran-like (1: x, 2: y) or not (0-based; 1: y, 2: x)
[2356]	1013	! vals: list of values correspondant to the coordinates
	1014	! dx, dy: shape of the 2D matrix
	1015	! NOval: Value to assign when there is no coordinate
	1016	! matrix: resultant matrix
	1017
	1018	fname = 'from_ptlist_2DRKmatrix'
	1019
	1020	matrix = NOval
	1021
[2359]	1022	IF (Flike) THEN
	1023	DO iv=1, Npts
	1024	i = pts(iv,1)
	1025	j = pts(iv,2)
	1026	matrix(i,j) = vals(iv)
	1027	END DO
	1028	ELSE
	1029	DO iv=1, Npts
	1030	i = pts(iv,2)+1
	1031	j = pts(iv,1)+1
	1032	matrix(i,j) = vals(iv)
	1033	END DO
	1034	END IF
[2356]	1035
	1036	RETURN
	1037
	1038	END SUBROUTINE from_ptlist_2DRKmatrix
	1039
[2359]	1040	SUBROUTINE from_ptlist_2DRKNmatrix(Npts, pts, Flike, Nvals, vals, dx, dy, NOval, Nmatrix)
[2357]	1041	! Subroutine to construct N 2D RK matrix from a list of values accompaigned by a list of grid-point
	1042	! coordinates
	1043
	1044	IMPLICIT NONE
	1045
	1046	INTEGER, INTENT(in) :: Npts, dx, dy, Nvals
[2359]	1047	LOGICAL, INTENT(in) :: Flike
[2357]	1048	REAL(r_k), INTENT(in) :: NOval
	1049	INTEGER, DIMENSION(Npts,2), INTENT(in) :: pts
	1050	REAL(r_k), DIMENSION(Npts,Nvals), INTENT(in) :: vals
	1051	REAL(r_k), DIMENSION(dx,dy,Nvals), INTENT(out) :: Nmatrix
	1052
	1053	! Local
	1054	INTEGER :: iv, i, j, ic
	1055
	1056	!!!!!!! Variables
	1057	! Npts: Number of values of the list
	1058	! pts: 2D matrix coordinates of the values
[2359]	1059	! Flike: whether input is Fortran-like (1: x, 2: y) or not (0-based; 1: y, 2: x)
[2357]	1060	! Nvals: number of different values
	1061	! vals: list of N-values correspondant to the coordinates
	1062	! dx, dy: shape of the 2D matrix
	1063	! NOval: Value to assign when there is no coordinate
	1064	! Nmatrix: resultant N-matrix
	1065
	1066	fname = 'from_ptlist_2DRKNmatrix'
	1067
	1068	Nmatrix = NOval
	1069
[2359]	1070	IF (Flike) THEN
	1071	DO iv=1, Npts
	1072	i = pts(iv,1)
	1073	j = pts(iv,2)
	1074	Nmatrix(i,j,:) = vals(iv,:)
	1075	END DO
	1076	ELSE
	1077	DO iv=1, Npts
	1078	i = pts(iv,2)+1
	1079	j = pts(iv,1)+1
	1080	Nmatrix(i,j,:) = vals(iv,:)
	1081	END DO
	1082	END IF
[2357]	1083
	1084	RETURN
	1085
	1086	END SUBROUTINE from_ptlist_2DRKNmatrix
	1087
[2359]	1088	SUBROUTINE from_coordlist_2DRKmatrix(Npts, coords, Flike, xcoord, ycoord, vals, dx, dy, NOval, &
	1089	matrix)
[2357]	1090	! Subroutine to construct a 2D RK matrix from a list of values accompaigned by a list of coordinates
	1091	! to find i,j grid-point coordinates by minimum distance
	1092
	1093	IMPLICIT NONE
	1094
	1095	INTEGER, INTENT(in) :: Npts, dx, dy
[2359]	1096	LOGICAL, INTENT(in) :: Flike
[2357]	1097	REAL(r_k), INTENT(in) :: NOval
	1098	REAL(r_k), DIMENSION(Npts,2), INTENT(in) :: coords
	1099	REAL(r_k), DIMENSION(Npts), INTENT(in) :: vals
	1100	REAL(r_k), DIMENSION(dx,dy), INTENT(in) :: xcoord, ycoord
	1101	REAL(r_k), DIMENSION(dx,dy), INTENT(out) :: matrix
	1102
	1103	! Local
	1104	INTEGER :: iv, i, j
	1105	REAL(r_k) :: xi, yi
	1106	REAL(r_k), DIMENSION(dx,dy) :: diff
	1107	INTEGER, DIMENSION(2) :: minpt
	1108
	1109	!!!!!!! Variables
	1110	! Npts: Number of values of the list
	1111	! coords: 2D coordinates of the values
[2359]	1112	! Flike: whether input is Fortran-like (1: x, 2: y) or not (0-based; 1: y, 2: x)
[2357]	1113	! vals: list of values correspondant to the coordinates
	1114	! xcoord, ycoord: matrix of values correspondant to the each coordinate
	1115	! dx, dy: shape of the 2D matrix
	1116	! NOval: Value to assign when there is no coordinate
	1117	! matrix: resultant matrix
	1118
	1119	fname = 'from_coordlist_2DRKmatrix'
	1120
	1121	matrix = NOval
	1122
[2359]	1123	IF (Flike) THEN
	1124	DO iv=1, Npts
	1125	xi = coords(iv,1)
	1126	yi = coords(iv,2)
	1127	diff = SQRT((xcoord-xi)2+(ycoord-yi)2)
	1128	minpt = MINLOC(diff)
	1129	matrix(minpt(1),minpt(2)) = vals(iv)
	1130	END DO
	1131	ELSE
	1132	DO iv=1, Npts
	1133	xi = coords(iv,2)
	1134	yi = coords(iv,1)
	1135	diff = SQRT((xcoord-xi)2+(ycoord-yi)2)
	1136	minpt = MINLOC(diff)
	1137	matrix(minpt(1),minpt(2)) = vals(iv)
	1138	END DO
	1139	END IF
[2357]	1140
	1141	RETURN
	1142
	1143	END SUBROUTINE from_coordlist_2DRKmatrix
	1144
[715]	1145	END MODULE module_generic

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format