Ignore:
Timestamp:
May 30, 2017 5:47:52 PM (4 years ago)
Author:
suehring
Message:

Adjustments according new topography and surface-modelling concept implemented

File:
1 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/SOURCE/wind_turbine_model_mod.f90

    r2153 r2232  
    2020! Current revisions:
    2121! -----------------
    22 !
     22! Adjustments to new topography concept
    2323!
    2424! Former revisions:
     
    8686    USE indices,                                                               &
    8787        ONLY:  nbgp, nx, nxl, nxlg, nxr, nxrg, ny, nyn, nyng, nys, nysg, nz,   &
    88                nzb, nzb_u_inner, nzb_v_inner, nzb_w_inner, nzt
     88               nzb, nzt, wall_flags_0
    8989
    9090    USE kinds
     
    13581358       INTEGER(iwp) ::  ii, jj, kk       !<
    13591359   
     1360       REAL(wp)     ::  flag               !< flag to mask topography grid points
    13601361       REAL(wp)     ::  sin_rot, cos_rot   !<
    13611362       REAL(wp)     ::  sin_yaw, cos_yaw   !<
     
    19711972                DO j = MAX( nys, j_hub(inot) - j_smear(inot) ),                &
    19721973                        MIN( nyn, j_hub(inot) + j_smear(inot) )
    1973                    DO k = MAX( nzb_u_inner(j,i)+1, k_hub(inot) - k_smear(inot) ), &
    1974                                 k_hub(inot) + k_smear(inot)
     1974!                    DO k = MAX( nzb_u_inner(j,i)+1, k_hub(inot) - k_smear(inot) ), &
     1975!                                 k_hub(inot) + k_smear(inot)
     1976                   DO  k = nzb+1, k_hub(inot) + k_smear(inot)
    19751977                      DO ring = 1, nrings(inot)
    19761978                         DO rseg = 1, nsegs(ring,inot)
     1979!
     1980!--                         Predetermine flag to mask topography
     1981                            flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_0(k,j,i), 1 ) )
     1982     
    19771983!
    19781984!--                         Determine the square of the distance between the
     
    19952001!--                         squared distances <= eps_min2:
    19962002                            IF ( dist_u_3d <= eps_min2 ) THEN
    1997                             thrust(k,j,i) = thrust(k,j,i) +                    &
    1998                                             thrust_ring(ring,rseg) *           &
    1999                                             ( ( pol_a * dist_u_3d - pol_b ) *  &
    2000                                              dist_u_3d + 1.0_wp ) * eps_factor
     2003                               thrust(k,j,i) = thrust(k,j,i) +                     &
     2004                                               thrust_ring(ring,rseg) *            &
     2005                                               ( ( pol_a * dist_u_3d - pol_b ) *   &
     2006                                                dist_u_3d + 1.0_wp ) * eps_factor *&
     2007                                                                       flag
    20012008                            ENDIF
    20022009                            IF ( dist_v_3d <= eps_min2 ) THEN
    2003                             torque_y(k,j,i) = torque_y(k,j,i) +                &
    2004                                               torque_ring_y(ring,rseg) *       &
    2005                                               ( ( pol_a * dist_v_3d - pol_b ) *&
    2006                                                dist_v_3d + 1.0_wp ) * eps_factor
     2010                               torque_y(k,j,i) = torque_y(k,j,i) +                    &
     2011                                                 torque_ring_y(ring,rseg) *           &
     2012                                                 ( ( pol_a * dist_v_3d - pol_b ) *    &
     2013                                                  dist_v_3d + 1.0_wp ) * eps_factor * &
     2014                                                                         flag
    20072015                            ENDIF
    20082016                            IF ( dist_w_3d <= eps_min2 ) THEN
    2009                             torque_z(k,j,i) = torque_z(k,j,i) +                &
    2010                                               torque_ring_z(ring,rseg) *       &
    2011                                               ( ( pol_a * dist_w_3d - pol_b ) *&
    2012                                                dist_w_3d + 1.0_wp ) * eps_factor
     2017                               torque_z(k,j,i) = torque_z(k,j,i) +                    &
     2018                                                 torque_ring_z(ring,rseg) *           &
     2019                                                 ( ( pol_a * dist_w_3d - pol_b ) *    &
     2020                                                  dist_w_3d + 1.0_wp ) * eps_factor * &
     2021                                                                         flag
    20132022                            ENDIF
    20142023
     
    20182027!
    20192028!--                   Rotation of force components:
    2020                       rot_tend_x(k,j,i) = rot_tend_x(k,j,i) +                  &
     2029                      rot_tend_x(k,j,i) = rot_tend_x(k,j,i) + (                &
    20212030                                      thrust(k,j,i)*rotx(inot,1) +             &
    20222031                                      torque_y(k,j,i)*roty(inot,1) +           &
    2023                                       torque_z(k,j,i)*rotz(inot,1)
     2032                                      torque_z(k,j,i)*rotz(inot,1)             &
     2033                                                              ) * flag
    20242034                               
    2025                       rot_tend_y(k,j,i) = rot_tend_y(k,j,i) +                  &
     2035                      rot_tend_y(k,j,i) = rot_tend_y(k,j,i) + (                &
    20262036                                      thrust(k,j,i)*rotx(inot,2) +             &
    20272037                                      torque_y(k,j,i)*roty(inot,2) +           &
    2028                                       torque_z(k,j,i)*rotz(inot,2)
     2038                                      torque_z(k,j,i)*rotz(inot,2)             &
     2039                                                              ) * flag
    20292040                               
    2030                       rot_tend_z(k,j,i) = rot_tend_z(k,j,i) +                  &
     2041                      rot_tend_z(k,j,i) = rot_tend_z(k,j,i) + (                &
    20312042                                      thrust(k,j,i)*rotx(inot,3) +             &
    20322043                                      torque_y(k,j,i)*roty(inot,3) +           &
    2033                                       torque_z(k,j,i)*rotz(inot,3)                               
     2044                                      torque_z(k,j,i)*rotz(inot,3)             &
     2045                                                              ) * flag                   
    20342046
    20352047                   ENDDO        ! End of loop over k
     
    24242436             DO  i = nxlg, nxrg
    24252437                DO  j = nysg, nyng
    2426                    DO  k = nzb_u_inner(j,i)+1, k_hub(1) + k_smear(1)
     2438                   DO  k = nzb+1, k_hub(1) + k_smear(1)
    24272439!
    24282440!--                   Calculate the thrust generated by the nacelle and the tower
     
    24322444                                         SIGN( u(k,j,i)**2 , u(k,j,i) )
    24332445                                                   
    2434                       tend(k,j,i) = tend(k,j,i) - rot_tend_x(k,j,i)            &
    2435                                   - tend_nac_x - tend_tow_x
     2446                      tend(k,j,i) = tend(k,j,i) + ( - rot_tend_x(k,j,i)        &
     2447                                  - tend_nac_x - tend_tow_x )                  &
     2448                                      * MERGE( 1.0_wp, 0.0_wp,                 &
     2449                                               BTEST( wall_flags_0(k,j,i), 1 ) )
    24362450                   ENDDO
    24372451                ENDDO
     
    24452459             DO  i = nxlg, nxrg
    24462460                DO  j = nysg, nyng
    2447                    DO  k = nzb_v_inner(j,i)+1, k_hub(1) + k_smear(1)
     2461                   DO  k = nzb+1, k_hub(1) + k_smear(1)
    24482462                      tend_nac_y = 0.5_wp * nac_cd_surf(k,j,i) *               &
    24492463                                         SIGN( v(k,j,i)**2 , v(k,j,i) )     
    24502464                      tend_tow_y   = 0.5_wp * tow_cd_surf(k,j,i) *             &
    24512465                                         SIGN( v(k,j,i)**2 , v(k,j,i) )                     
    2452                       tend(k,j,i) = tend(k,j,i) - rot_tend_y(k,j,i)            &
    2453                                   - tend_nac_y - tend_tow_y
     2466                      tend(k,j,i) = tend(k,j,i) + ( - rot_tend_y(k,j,i)        &
     2467                                  - tend_nac_y - tend_tow_y )                  &
     2468                                      * MERGE( 1.0_wp, 0.0_wp,                 &
     2469                                               BTEST( wall_flags_0(k,j,i), 2 ) )
    24542470                   ENDDO
    24552471                ENDDO
     
    24632479             DO  i = nxlg, nxrg
    24642480                DO  j = nysg, nyng
    2465                    DO  k = nzb_w_inner(j,i)+1,  k_hub(1) + k_smear(1)
    2466                       tend(k,j,i) = tend(k,j,i) - rot_tend_z(k,j,i)
     2481                   DO  k = nzb+1,  k_hub(1) + k_smear(1)
     2482                      tend(k,j,i) = tend(k,j,i) - rot_tend_z(k,j,i)            &
     2483                                      * MERGE( 1.0_wp, 0.0_wp,                 &
     2484                                               BTEST( wall_flags_0(k,j,i), 3 ) )
    24672485                   ENDDO
    24682486                ENDDO
     
    25062524          IF ( simulated_time >= time_turbine_on )  THEN
    25072525
    2508              DO  k = nzb_u_inner(j,i)+1,  k_hub(1) + k_smear(1)
     2526             DO  k = nzb+1,  k_hub(1) + k_smear(1)
    25092527!
    25102528!--             Calculate the thrust generated by the nacelle and the tower
     
    25132531                tend_tow_x   = 0.5_wp * tow_cd_surf(k,j,i) *                   &
    25142532                                   SIGN( u(k,j,i)**2 , u(k,j,i) )
    2515                 tend(k,j,i) = tend(k,j,i) - rot_tend_x(k,j,i)                  &
    2516                             - tend_nac_x - tend_tow_x
     2533                tend(k,j,i) = tend(k,j,i) + ( - rot_tend_x(k,j,i)              &
     2534                            - tend_nac_x - tend_tow_x )                        &
     2535                                      * MERGE( 1.0_wp, 0.0_wp,                 &
     2536                                               BTEST( wall_flags_0(k,j,i), 1 ) )
    25172537             ENDDO
    25182538          ENDIF
     
    25222542!--       Apply the y-component of the force to the v-component of the flow:
    25232543          IF ( simulated_time >= time_turbine_on )  THEN
    2524              DO  k = nzb_v_inner(j,i)+1,  k_hub(1) + k_smear(1)
     2544             DO  k = nzb+1,  k_hub(1) + k_smear(1)
    25252545                tend_nac_y = 0.5_wp * nac_cd_surf(k,j,i) *                     &
    25262546                                   SIGN( v(k,j,i)**2 , v(k,j,i) )     
    25272547                tend_tow_y   = 0.5_wp * tow_cd_surf(k,j,i) *                   &
    25282548                                   SIGN( v(k,j,i)**2 , v(k,j,i) )                     
    2529                 tend(k,j,i) = tend(k,j,i) - rot_tend_y(k,j,i)                  &
    2530                             - tend_nac_y - tend_tow_y
     2549                tend(k,j,i) = tend(k,j,i) + ( - rot_tend_y(k,j,i)              &
     2550                            - tend_nac_y - tend_tow_y )                        &
     2551                                      * MERGE( 1.0_wp, 0.0_wp,                 &
     2552                                               BTEST( wall_flags_0(k,j,i), 2 ) )
    25312553             ENDDO
    25322554          ENDIF
     
    25362558!--       Apply the z-component of the force to the w-component of the flow:
    25372559          IF ( simulated_time >= time_turbine_on )  THEN
    2538              DO  k = nzb_w_inner(j,i)+1,  k_hub(1) + k_smear(1)
    2539                 tend(k,j,i) = tend(k,j,i) - rot_tend_z(k,j,i)
     2560             DO  k = nzb+1,  k_hub(1) + k_smear(1)
     2561                tend(k,j,i) = tend(k,j,i) - rot_tend_z(k,j,i)                  &
     2562                                      * MERGE( 1.0_wp, 0.0_wp,                 &
     2563                                               BTEST( wall_flags_0(k,j,i), 3 ) )
    25402564             ENDDO
    25412565          ENDIF
Note: See TracChangeset for help on using the changeset viewer.