[1] | 1 | SUBROUTINE asselin_filter |
---|
| 2 | |
---|
[3] | 3 | !------------------------------------------------------------------------------! |
---|
[1] | 4 | ! Actual revisions: |
---|
| 5 | ! ----------------- |
---|
| 6 | ! |
---|
| 7 | ! |
---|
| 8 | ! Former revisions: |
---|
| 9 | ! --------------------- |
---|
[3] | 10 | ! $Id: asselin_filter.f90 4 2007-02-13 11:33:16Z raasch $ |
---|
| 11 | ! RCS Log replace by Id keyword, revision history cleaned up |
---|
| 12 | ! |
---|
[1] | 13 | ! Revision 1.8 2004/01/30 10:14:02 raasch |
---|
| 14 | ! Scalar lower k index nzb replaced by 2d-array nzb_2d |
---|
| 15 | ! |
---|
| 16 | ! Revision 1.1 2002/05/02 13:43:53 13:43:53 raasch (Siegfried Raasch) |
---|
| 17 | ! Initial revision |
---|
| 18 | ! |
---|
| 19 | ! |
---|
| 20 | ! Description: |
---|
| 21 | ! ------------- |
---|
| 22 | ! Time filter needed for the leap-frog method |
---|
[3] | 23 | !------------------------------------------------------------------------------! |
---|
[1] | 24 | |
---|
| 25 | USE arrays_3d |
---|
| 26 | USE control_parameters |
---|
| 27 | USE cpulog |
---|
| 28 | USE indices |
---|
| 29 | USE interfaces |
---|
| 30 | |
---|
| 31 | IMPLICIT NONE |
---|
| 32 | |
---|
| 33 | INTEGER :: i, j, k |
---|
| 34 | |
---|
| 35 | |
---|
| 36 | CALL cpu_log( log_point(9), 'timefilter', 'start' ) |
---|
| 37 | |
---|
| 38 | ! |
---|
| 39 | !-- Return to the calling routine, if time filter is not to be applied |
---|
| 40 | IF ( asselin_filter_factor == 0.0 ) RETURN |
---|
| 41 | |
---|
| 42 | ! |
---|
| 43 | !-- Apply the time filter |
---|
| 44 | #if defined( __ibm ) |
---|
| 45 | !$OMP PARALLEL PRIVATE (i,j,k) |
---|
| 46 | !$OMP DO |
---|
| 47 | DO i = nxl-1, nxr+1 |
---|
| 48 | DO j = nys-1, nyn+1 |
---|
| 49 | |
---|
| 50 | DO k = nzb_2d(j,i), nzt+1 |
---|
| 51 | u(k,j,i) = u(k,j,i) + asselin_filter_factor * & |
---|
| 52 | ( u_p(k,j,i) - 2.0 * u(k,j,i) + u_m(k,j,i) ) |
---|
| 53 | v(k,j,i) = v(k,j,i) + asselin_filter_factor * & |
---|
| 54 | ( v_p(k,j,i) - 2.0 * v(k,j,i) + v_m(k,j,i) ) |
---|
| 55 | w(k,j,i) = w(k,j,i) + asselin_filter_factor * & |
---|
| 56 | ( w_p(k,j,i) - 2.0 * w(k,j,i) + w_m(k,j,i) ) |
---|
| 57 | ENDDO |
---|
| 58 | |
---|
| 59 | IF ( scalar_advec /= 'bc-scheme' ) THEN |
---|
| 60 | DO k = nzb_2d(j,i), nzt+1 |
---|
| 61 | pt(k,j,i) = pt(k,j,i) + asselin_filter_factor * & |
---|
| 62 | ( pt_p(k,j,i) - 2.0 * pt(k,j,i) + pt_m(k,j,i) ) |
---|
| 63 | ENDDO |
---|
| 64 | ENDIF |
---|
| 65 | |
---|
| 66 | IF ( .NOT. constant_diffusion .AND. scalar_advec /= 'bc-scheme' ) & |
---|
| 67 | THEN |
---|
| 68 | DO k = nzb_2d(j,i), nzt+1 |
---|
| 69 | e(k,j,i) = e(k,j,i) + asselin_filter_factor * & |
---|
| 70 | ( e_p(k,j,i) - 2.0 * e(k,j,i) + e_m(k,j,i) ) |
---|
| 71 | ENDDO |
---|
| 72 | ENDIF |
---|
| 73 | |
---|
| 74 | IF ( ( moisture .OR. passive_scalar ) .AND. & |
---|
| 75 | scalar_advec /= 'bc-scheme' ) THEN |
---|
| 76 | DO k = nzb_2d(j,i), nzt+1 |
---|
| 77 | q(k,j,i) = q(k,j,i) + asselin_filter_factor * & |
---|
| 78 | ( q_p(k,j,i) - 2.0 * q(k,j,i) + q_m(k,j,i) ) |
---|
| 79 | ENDDO |
---|
| 80 | ENDIF |
---|
| 81 | |
---|
| 82 | ENDDO |
---|
| 83 | ENDDO |
---|
| 84 | !$OMP END PARALLEL |
---|
| 85 | #else |
---|
| 86 | u = u + asselin_filter_factor * ( u_p - 2.0 * u + u_m ) |
---|
| 87 | v = v + asselin_filter_factor * ( v_p - 2.0 * v + v_m ) |
---|
| 88 | w = w + asselin_filter_factor * ( w_p - 2.0 * w + w_m ) |
---|
| 89 | |
---|
| 90 | IF ( scalar_advec /= 'bc-scheme' ) THEN |
---|
| 91 | pt = pt + asselin_filter_factor * ( pt_p - 2.0 * pt + pt_m ) |
---|
| 92 | ENDIF |
---|
| 93 | |
---|
| 94 | IF ( .NOT. constant_diffusion .AND. scalar_advec /= 'bc-scheme' ) THEN |
---|
| 95 | e = e + asselin_filter_factor * ( e_p - 2.0 * e + e_m ) |
---|
| 96 | ENDIF |
---|
| 97 | |
---|
| 98 | IF ( ( moisture .OR. passive_scalar ) .AND. & |
---|
| 99 | scalar_advec /= 'bc-scheme' ) THEN |
---|
| 100 | q = q + asselin_filter_factor * ( q_p - 2.0 * q + q_m ) |
---|
| 101 | ENDIF |
---|
| 102 | #endif |
---|
| 103 | |
---|
| 104 | CALL cpu_log( log_point(9), 'timefilter', 'stop' ) |
---|
| 105 | |
---|
| 106 | END SUBROUTINE asselin_filter |
---|