Parent Directory | Revision Log
Clarified the logic in filtreg by creating a procedure filtreg_hemisph. It was terrible with a loop on hemispheres and tests on hemisphere inside the loop, plus maddening indirections on latitude bounds, plus repeated code. Went from 126 lines to much clearer 74 + 32 = 106 lines. In module inifilr_m, finally made the arrays matrice[uv][ns], matrinv[ns] dynamic (following LMDZ). Changed the lower bound of matriceun and matrinvn in the 3rd dimension: 2 instead of 1, the index 1 was not defined (nor used). In module inifilr_m, changed the bounds of matriceus and matrinvs in the 3rd dimension: jfiltsu:jjm instead of 1:jjm - jfiltsu + 1. Changed the bounds of matricevs in the 3rd dimension: jfiltsv:jjm instead of 1:jjm - jfiltsv + 1. It is a little simpler and clearer this way in procedure inifilr.
1 | module filtreg_hemisph_m |
2 | |
3 | implicit none |
4 | |
5 | contains |
6 | |
7 | subroutine filtreg_hemisph(champ, sdd1, sdd2, matri) |
8 | |
9 | USE dimens_m, ONLY: iim |
10 | |
11 | REAL, intent(inout):: champ(:, :, :) ! (iim + 1, :, :) |
12 | REAL, intent(in):: sdd1(:), sdd2(:) ! (iim) |
13 | real, intent(in), dimension(:, :, :):: matri ! (iim, iim, :) |
14 | |
15 | ! Local: |
16 | integer l, j |
17 | |
18 | !----------------------------------------------------------------- |
19 | |
20 | DO l = 1, size(champ, 3) |
21 | DO j = 1, size(champ, 2) |
22 | champ(:iim, j, l) = champ(:iim, j, l) * sdd1 |
23 | champ(:iim, j, l) = (champ(:iim, j, l) & |
24 | + matmul(matri(:, :, j), champ(:iim, j, l))) * sdd2 |
25 | champ(iim + 1, j, l) = champ(1, j, l) |
26 | END DO |
27 | END DO |
28 | |
29 | end subroutine filtreg_hemisph |
30 | |
31 | end module filtreg_hemisph_m |
ViewVC Help | |
Powered by ViewVC 1.1.21 |