Changeset 2553 for palm/trunk/SOURCE/wind_turbine_model_mod.f90
 Timestamp:
 Oct 18, 2017 8:03:45 AM (4 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

palm/trunk/SOURCE/wind_turbine_model_mod.f90
r2410 r2553 26 26 !  27 27 ! $Id$ 28 ! Bugfix of vertical loop in wtm_tendencies to account for different turbine 29 ! heights, bugfix of the interpolation of the ucomponent concerning the 30 ! vertical index and further small adjustments of the programming style 31 ! 32 ! 2410 20170906 08:16:30Z Giersch 28 33 ! Revise error message PA0462 29 34 ! … … 93 98 94 99 USE arrays_3d, & 95 ONLY: dd2zu, tend, u, v, w, zu, zw 96 97 USE constants 100 ONLY: tend, u, v, w, zu, zw 101 102 USE constants, & 103 ONLY: pi 98 104 99 105 USE control_parameters, & … … 140 146 141 147 REAL(wp), DIMENSION(1:100) :: dtow = 0.0_wp !< tower diameter [m] 142 REAL(wp), DIMENSION(1:100) :: omega_rot = 0. 0_wp !< inital or constant rotor speed [rad/s]148 REAL(wp), DIMENSION(1:100) :: omega_rot = 0.9_wp !< inital or constant rotor speed [rad/s] 143 149 REAL(wp), DIMENSION(1:100) :: phi_yaw = 0.0_wp !< yaw angle [degree] ( clockwise, 0 = facing west ) 144 150 REAL(wp), DIMENSION(1:100) :: pitch_add = 0.0_wp !< constant pitch angle … … 240 246 241 247 REAL(wp) :: cur_r !< 242 REAL(wp) :: delta_t !< tangential segment length243 248 REAL(wp) :: phi_rotor !< 244 249 REAL(wp) :: pre_factor !< … … 383 388 INTERFACE wtm_forces 384 389 MODULE PROCEDURE wtm_forces 390 END INTERFACE wtm_forces 391 392 INTERFACE wtm_yawcontrol 385 393 MODULE PROCEDURE wtm_yawcontrol 386 END INTERFACE wtm_ forces394 END INTERFACE wtm_yawcontrol 387 395 388 396 INTERFACE wtm_rotate_rotor … … 402 410 403 411 PUBLIC wtm_check_parameters, wtm_forces, wtm_init, wtm_init_arrays, & 404 wtm_parin, wtm_tendencies, w tm_tendencies_ij, wind_turbine412 wtm_parin, wtm_tendencies, wind_turbine 405 413 406 414 CONTAINS … … 524 532 delta_t_factor = segment_length 525 533 delta_r_init = delta_r_factor * MIN( dx, dy, dz) 526 delta_t = delta_t_factor * MIN( dx, dy, dz)527 534 528 535 DO inot = 1, nturbines … … 1484 1491 ii = rbx(ring,rseg) * ddx 1485 1492 jj = ( rby(ring,rseg)  0.5_wp * dy ) * ddy 1486 kk = ( rbz(ring,rseg) 0.5_wp * dz ) / dz1493 kk = ( rbz(ring,rseg) + 0.5_wp * dz ) / dz 1487 1494 ! 1488 1495 ! Interpolate only if all required information is available on … … 1710 1717 ! 1711 1718 ! In the following the 3Dvelocity field is projected its 1712 ! components perpe dicular and parallel to the rotor area1719 ! components perpendicular and parallel to the rotor area 1713 1720 ! The calculation of forces will be done in the rotor 1714 1721 ! coordinates y' and z. … … 2230 2237 SUBROUTINE wtm_yawcontrol( inot ) 2231 2238 2232 USE constants2233 2239 USE kinds 2234 2240 … … 2461 2467 DO i = nxlg, nxrg 2462 2468 DO j = nysg, nyng 2463 DO k = nzb+1, k_hub(1) + k_smear(1)2469 DO k = nzb+1, MAXVAL(k_hub) + MAXVAL(k_smear) 2464 2470 ! 2465 2471 ! Calculate the thrust generated by the nacelle and the tower 2466 2472 tend_nac_x = 0.5_wp * nac_cd_surf(k,j,i) * & 2467 2473 SIGN( u(k,j,i)**2 , u(k,j,i) ) 2468 2474 tend_tow_x = 0.5_wp * tow_cd_surf(k,j,i) * & 2469 2475 SIGN( u(k,j,i)**2 , u(k,j,i) ) 2470 2476 2471 2477 tend(k,j,i) = tend(k,j,i) + (  rot_tend_x(k,j,i) & 2472 2478  tend_nac_x  tend_tow_x ) & … … 2484 2490 DO i = nxlg, nxrg 2485 2491 DO j = nysg, nyng 2486 DO k = nzb+1, k_hub(1) + k_smear(1)2492 DO k = nzb+1, MAXVAL(k_hub) + MAXVAL(k_smear) 2487 2493 tend_nac_y = 0.5_wp * nac_cd_surf(k,j,i) * & 2488 2494 SIGN( v(k,j,i)**2 , v(k,j,i) ) … … 2504 2510 DO i = nxlg, nxrg 2505 2511 DO j = nysg, nyng 2506 DO k = nzb+1, k_hub(1) + k_smear(1)2512 DO k = nzb+1, MAXVAL(k_hub) + MAXVAL(k_smear) 2507 2513 tend(k,j,i) = tend(k,j,i)  rot_tend_z(k,j,i) & 2508 2514 * MERGE( 1.0_wp, 0.0_wp, & … … 2548 2554 ! Apply the xcomponent of the force to the ucomponent of the flow: 2549 2555 IF ( simulated_time >= time_turbine_on ) THEN 2550 2551 DO k = nzb+1, k_hub(1) + k_smear(1) 2556 DO k = nzb+1, MAXVAL(k_hub) + MAXVAL(k_smear) 2552 2557 ! 2553 2558 ! Calculate the thrust generated by the nacelle and the tower … … 2559 2564  tend_nac_x  tend_tow_x ) & 2560 2565 * MERGE( 1.0_wp, 0.0_wp, & 2561 2566 BTEST( wall_flags_0(k,j,i), 1 ) ) 2562 2567 ENDDO 2563 2568 ENDIF … … 2567 2572 ! Apply the ycomponent of the force to the vcomponent of the flow: 2568 2573 IF ( simulated_time >= time_turbine_on ) THEN 2569 DO k = nzb+1, k_hub(1) + k_smear(1)2574 DO k = nzb+1, MAXVAL(k_hub) + MAXVAL(k_smear) 2570 2575 tend_nac_y = 0.5_wp * nac_cd_surf(k,j,i) * & 2571 2576 SIGN( v(k,j,i)**2 , v(k,j,i) ) … … 2576 2581 * MERGE( 1.0_wp, 0.0_wp, & 2577 2582 BTEST( wall_flags_0(k,j,i), 2 ) ) 2578 ENDDO2583 ENDDO 2579 2584 ENDIF 2580 2585 … … 2583 2588 ! Apply the zcomponent of the force to the wcomponent of the flow: 2584 2589 IF ( simulated_time >= time_turbine_on ) THEN 2585 DO k = nzb+1, k_hub(1) + k_smear(1)2590 DO k = nzb+1, MAXVAL(k_hub) + MAXVAL(k_smear) 2586 2591 tend(k,j,i) = tend(k,j,i)  rot_tend_z(k,j,i) & 2587 2592 * MERGE( 1.0_wp, 0.0_wp, &
Note: See TracChangeset
for help on using the changeset viewer.