Changeset 3712 for palm/trunk/SOURCE
- Timestamp:
- Feb 1, 2019 3:37:28 PM (4 years ago)
- File:
-
- 1 edited
-
palm/trunk/SOURCE/urban_surface_mod.f90 (modified) (248 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/urban_surface_mod.f90
r3710 r3712 15 15 ! PALM. If not, see <http://www.gnu.org/licenses/>. 16 16 ! 17 ! Copyright 2015-201 8Czech Technical University in Prague18 ! Copyright 2015-201 8Institute of Computer Science of the17 ! Copyright 2015-2019 Czech Technical University in Prague 18 ! Copyright 2015-2019 Institute of Computer Science of the 19 19 ! Czech Academy of Sciences, Prague 20 20 ! Copyright 1997-2019 Leibniz Universitaet Hannover … … 28 28 ! ----------------- 29 29 ! $Id$ 30 ! Formatting and clean-up (rvtils) 31 ! 32 ! 3710 2019-01-30 18:11:19Z suehring 30 33 ! Check if building type is set within a valid range. 31 34 ! … … 477 480 478 481 REAL(wp), PARAMETER :: & 479 b_ch = 6.04_wp, & !Clapp & Hornberger exponent480 lambda_h_green_dry = 0.19_wp, & !heat conductivity for dry soil481 lambda_h_green_sm = 3.44_wp, & !heat conductivity of the soil matrix482 lambda_h_water = 0.57_wp, & !heat conductivity of water483 psi_sat = -0.388_wp, & !soil matrix potential at saturation484 rho_c_soil = 2.19E6_wp, & !volumetric heat capacity of soil485 rho_c_water = 4.20E6_wp !, & !volumetric heat capacity of water482 b_ch = 6.04_wp, & !< Clapp & Hornberger exponent 483 lambda_h_green_dry = 0.19_wp, & !< heat conductivity for dry soil 484 lambda_h_green_sm = 3.44_wp, & !< heat conductivity of the soil matrix 485 lambda_h_water = 0.57_wp, & !< heat conductivity of water 486 psi_sat = -0.388_wp, & !< soil matrix potential at saturation 487 rho_c_soil = 2.19E6_wp, & !< volumetric heat capacity of soil 488 rho_c_water = 4.20E6_wp !< volumetric heat capacity of water 486 489 ! m_max_depth = 0.0002_wp ! Maximum capacity of the water reservoir (m) 487 490 … … 489 492 !-- Soil parameters I alpha_vg, l_vg_green, n_vg, gamma_w_green_sat 490 493 REAL(wp), DIMENSION(0:3,1:7), PARAMETER :: soil_pars = RESHAPE( (/ & 491 3.83_wp, 1.250_wp, 1.38_wp, 6.94E-6_wp, & !1492 3.14_wp, -2.342_wp, 1.28_wp, 1.16E-6_wp, & !2493 0.83_wp, -0.588_wp, 1.25_wp, 0.26E-6_wp, & !3494 3.67_wp, -1.977_wp, 1.10_wp, 2.87E-6_wp, & !4495 2.65_wp, 2.500_wp, 1.10_wp, 1.74E-6_wp, & !5496 1.30_wp, 0.400_wp, 1.20_wp, 0.93E-6_wp, & !6497 0.00_wp, 0.00_wp, 0.00_wp, 0.57E-6_wp & !7494 3.83_wp, 1.250_wp, 1.38_wp, 6.94E-6_wp, & !< soil 1 495 3.14_wp, -2.342_wp, 1.28_wp, 1.16E-6_wp, & !< soil 2 496 0.83_wp, -0.588_wp, 1.25_wp, 0.26E-6_wp, & !< soil 3 497 3.67_wp, -1.977_wp, 1.10_wp, 2.87E-6_wp, & !< soil 4 498 2.65_wp, 2.500_wp, 1.10_wp, 1.74E-6_wp, & !< soil 5 499 1.30_wp, 0.400_wp, 1.20_wp, 0.93E-6_wp, & !< soil 6 500 0.00_wp, 0.00_wp, 0.00_wp, 0.57E-6_wp & !< soil 7 498 501 /), (/ 4, 7 /) ) 499 502 … … 501 504 !-- Soil parameters II swc_sat, fc, wilt, swc_res 502 505 REAL(wp), DIMENSION(0:3,1:7), PARAMETER :: m_soil_pars = RESHAPE( (/ & 503 0.403_wp, 0.244_wp, 0.059_wp, 0.025_wp, & !1504 0.439_wp, 0.347_wp, 0.151_wp, 0.010_wp, & !2505 0.430_wp, 0.383_wp, 0.133_wp, 0.010_wp, & !3506 0.520_wp, 0.448_wp, 0.279_wp, 0.010_wp, & !4507 0.614_wp, 0.541_wp, 0.335_wp, 0.010_wp, & !5508 0.766_wp, 0.663_wp, 0.267_wp, 0.010_wp, & !6509 0.472_wp, 0.323_wp, 0.171_wp, 0.000_wp & !7506 0.403_wp, 0.244_wp, 0.059_wp, 0.025_wp, & !< soil 1 507 0.439_wp, 0.347_wp, 0.151_wp, 0.010_wp, & !< soil 2 508 0.430_wp, 0.383_wp, 0.133_wp, 0.010_wp, & !< soil 3 509 0.520_wp, 0.448_wp, 0.279_wp, 0.010_wp, & !< soil 4 510 0.614_wp, 0.541_wp, 0.335_wp, 0.010_wp, & !< soil 5 511 0.766_wp, 0.663_wp, 0.267_wp, 0.010_wp, & !< soil 6 512 0.472_wp, 0.323_wp, 0.171_wp, 0.000_wp & !< soil 7 510 513 /), (/ 4, 7 /) ) 511 512 !value 9999999.9_wp -> generic available or user-defined value must be set513 ! 514 REAL(wp) :: alpha_vangenuchten = 9999999.9_wp, & !< NAMELIST alpha_vg515 field_capacity = 9999999.9_wp, & !< NAMELIST fc516 hydraulic_conductivity = 9999999.9_wp, & !< NAMELIST gamma_w_green_sat517 lambda_h_green_sat = 0.0_wp, & !< heat conductivity for saturated soil518 l_vangenuchten = 9999999.9_wp, & !< NAMELIST l_vg519 n_vangenuchten = 9999999.9_wp, & !< NAMELIST n_vg520 residual_moisture = 9999999.9_wp, & !< NAMELIST m_res521 saturation_moisture = 9999999.9_wp, & !< NAMELIST m_sat522 wilting_point = 9999999.9_wp !, &!< NAMELIST m_wilt514 ! 515 !-- value 9999999.9_wp -> generic available or user-defined value must be set 516 !-- otherwise -> no generic variable and user setting is optional 517 REAL(wp) :: alpha_vangenuchten = 9999999.9_wp, & !< NAMELIST alpha_vg 518 field_capacity = 9999999.9_wp, & !< NAMELIST fc 519 hydraulic_conductivity = 9999999.9_wp, & !< NAMELIST gamma_w_green_sat 520 lambda_h_green_sat = 0.0_wp, & !< heat conductivity for saturated soil 521 l_vangenuchten = 9999999.9_wp, & !< NAMELIST l_vg 522 n_vangenuchten = 9999999.9_wp, & !< NAMELIST n_vg 523 residual_moisture = 9999999.9_wp, & !< NAMELIST m_res 524 saturation_moisture = 9999999.9_wp, & !< NAMELIST m_sat 525 wilting_point = 9999999.9_wp !< NAMELIST m_wilt 523 526 524 527 ! 525 528 !-- configuration parameters (they can be setup in PALM config) 526 529 LOGICAL :: usm_material_model = .TRUE. !< flag parameter indicating wheather the model of heat in materials is used 527 LOGICAL :: usm_anthropogenic_heat = .FALSE. !< flag parameter indicating wheather the anthropogenic heat sources (e.g.transportation) are used 530 LOGICAL :: usm_anthropogenic_heat = .FALSE. !< flag parameter indicating wheather the anthropogenic heat sources 531 !< (e.g.transportation) are used 528 532 LOGICAL :: force_radiation_call_l = .FALSE. !< flag parameter for unscheduled radiation model calls 529 533 LOGICAL :: indoor_model = .FALSE. !< whether to use the indoor model … … 540 544 ! 541 545 !-- Indices of input attributes for (above) ground floor level 542 INTEGER(iwp) :: ind_alb_wall_agfl = 65 !< index in input list for albedo_type of wall above ground floor level 543 INTEGER(iwp) :: ind_alb_wall_gfl = 32 !< index in input list for albedo_type of wall ground floor level 544 INTEGER(iwp) :: ind_alb_wall_r = 96 !< index in input list for albedo_type of wall roof 545 INTEGER(iwp) :: ind_alb_green_agfl = 83 !< index in input list for albedo_type of green above ground floor level 546 INTEGER(iwp) :: ind_alb_green_gfl = 50 !< index in input list for albedo_type of green ground floor level 547 INTEGER(iwp) :: ind_alb_green_r = 115 !< index in input list for albedo_type of green roof 548 INTEGER(iwp) :: ind_alb_win_agfl = 79 !< index in input list for albedo_type of window fraction above ground floor level 549 INTEGER(iwp) :: ind_alb_win_gfl = 46 !< index in input list for albedo_type of window fraction ground floor level 550 INTEGER(iwp) :: ind_alb_win_r = 110 !< index in input list for albedo_type of window fraction roof 551 INTEGER(iwp) :: ind_emis_wall_agfl = 64 !< index in input list for wall emissivity, above ground floor level 552 INTEGER(iwp) :: ind_emis_wall_gfl = 31 !< index in input list for wall emissivity, ground floor level 553 INTEGER(iwp) :: ind_emis_wall_r = 95 !< index in input list for wall emissivity, roof 554 INTEGER(iwp) :: ind_emis_green_agfl = 82 !< index in input list for green emissivity, above ground floor level 555 INTEGER(iwp) :: ind_emis_green_gfl = 49 !< index in input list for green emissivity, ground floor level 556 INTEGER(iwp) :: ind_emis_green_r = 114 !< index in input list for green emissivity, roof 557 INTEGER(iwp) :: ind_emis_win_agfl = 77 !< index in input list for window emissivity, above ground floor level 558 INTEGER(iwp) :: ind_emis_win_gfl = 44 !< index in input list for window emissivity, ground floor level 559 INTEGER(iwp) :: ind_emis_win_r = 108 !< index in input list for window emissivity, roof 560 INTEGER(iwp) :: ind_green_frac_w_agfl = 80 !< index in input list for green fraction on wall, above ground floor level 561 INTEGER(iwp) :: ind_green_frac_w_gfl = 47 !< index in input list for green fraction on wall, ground floor level 562 INTEGER(iwp) :: ind_green_frac_r_agfl = 112 !< index in input list for green fraction on roof, above ground floor level 563 INTEGER(iwp) :: ind_green_frac_r_gfl = 111 !< index in input list for green fraction on roof, ground floor level 564 INTEGER(iwp) :: ind_hc1_agfl = 58 !< index in input list for heat capacity at first wall layer, above ground floor level 565 INTEGER(iwp) :: ind_hc1_gfl = 25 !< index in input list for heat capacity at first wall layer, ground floor level 566 INTEGER(iwp) :: ind_hc1_wall_r = 89 !< index in input list for heat capacity at first wall layer, roof 567 INTEGER(iwp) :: ind_hc1_win_agfl = 71 !< index in input list for heat capacity at first window layer, above ground floor level 568 INTEGER(iwp) :: ind_hc1_win_gfl = 38 !< index in input list for heat capacity at first window layer, ground floor level 569 INTEGER(iwp) :: ind_hc1_win_r = 102 !< index in input list for heat capacity at first window layer, roof 570 INTEGER(iwp) :: ind_hc2_agfl = 59 !< index in input list for heat capacity at second wall layer, above ground floor level 571 INTEGER(iwp) :: ind_hc2_gfl = 26 !< index in input list for heat capacity at second wall layer, ground floor level 572 INTEGER(iwp) :: ind_hc2_wall_r = 90 !< index in input list for heat capacity at second wall layer, roof 573 INTEGER(iwp) :: ind_hc2_win_agfl = 72 !< index in input list for heat capacity at second window layer, above ground floor level 574 INTEGER(iwp) :: ind_hc2_win_gfl = 39 !< index in input list for heat capacity at second window layer, ground floor level 575 INTEGER(iwp) :: ind_hc2_win_r = 103 !< index in input list for heat capacity at second window layer, roof 576 INTEGER(iwp) :: ind_hc3_agfl = 60 !< index in input list for heat capacity at third wall layer, above ground floor level 577 INTEGER(iwp) :: ind_hc3_gfl = 27 !< index in input list for heat capacity at third wall layer, ground floor level 578 INTEGER(iwp) :: ind_hc3_wall_r = 91 !< index in input list for heat capacity at third wall layer, roof 579 INTEGER(iwp) :: ind_hc3_win_agfl = 73 !< index in input list for heat capacity at third window layer, above ground floor level 580 INTEGER(iwp) :: ind_hc3_win_gfl = 40 !< index in input list for heat capacity at third window layer, ground floor level 581 INTEGER(iwp) :: ind_hc3_win_r = 104 !< index in input list for heat capacity at third window layer, roof 582 INTEGER(iwp) :: ind_gflh = 17 !< index in input list for ground floor level height 583 INTEGER(iwp) :: ind_lai_r_agfl = 113 !< index in input list for LAI on roof, above ground floor level 546 INTEGER(iwp) :: ind_alb_wall_agfl = 65 !< index in input list for albedo_type of wall above ground floor level 547 INTEGER(iwp) :: ind_alb_wall_gfl = 32 !< index in input list for albedo_type of wall ground floor level 548 INTEGER(iwp) :: ind_alb_wall_r = 96 !< index in input list for albedo_type of wall roof 549 INTEGER(iwp) :: ind_alb_green_agfl = 83 !< index in input list for albedo_type of green above ground floor level 550 INTEGER(iwp) :: ind_alb_green_gfl = 50 !< index in input list for albedo_type of green ground floor level 551 INTEGER(iwp) :: ind_alb_green_r = 115 !< index in input list for albedo_type of green roof 552 INTEGER(iwp) :: ind_alb_win_agfl = 79 !< index in input list for albedo_type of window fraction 553 !< above ground floor level 554 INTEGER(iwp) :: ind_alb_win_gfl = 46 !< index in input list for albedo_type of window fraction ground floor level 555 INTEGER(iwp) :: ind_alb_win_r = 110 !< index in input list for albedo_type of window fraction roof 556 INTEGER(iwp) :: ind_emis_wall_agfl = 64 !< index in input list for wall emissivity, above ground floor level 557 INTEGER(iwp) :: ind_emis_wall_gfl = 31 !< index in input list for wall emissivity, ground floor level 558 INTEGER(iwp) :: ind_emis_wall_r = 95 !< index in input list for wall emissivity, roof 559 INTEGER(iwp) :: ind_emis_green_agfl = 82 !< index in input list for green emissivity, above ground floor level 560 INTEGER(iwp) :: ind_emis_green_gfl = 49 !< index in input list for green emissivity, ground floor level 561 INTEGER(iwp) :: ind_emis_green_r = 114 !< index in input list for green emissivity, roof 562 INTEGER(iwp) :: ind_emis_win_agfl = 77 !< index in input list for window emissivity, above ground floor level 563 INTEGER(iwp) :: ind_emis_win_gfl = 44 !< index in input list for window emissivity, ground floor level 564 INTEGER(iwp) :: ind_emis_win_r = 108 !< index in input list for window emissivity, roof 565 INTEGER(iwp) :: ind_green_frac_w_agfl = 80 !< index in input list for green fraction on wall, above ground floor level 566 INTEGER(iwp) :: ind_green_frac_w_gfl = 47 !< index in input list for green fraction on wall, ground floor level 567 INTEGER(iwp) :: ind_green_frac_r_agfl = 112 !< index in input list for green fraction on roof, above ground floor level 568 INTEGER(iwp) :: ind_green_frac_r_gfl = 111 !< index in input list for green fraction on roof, ground floor level 569 INTEGER(iwp) :: ind_hc1_agfl = 58 !< index in input list for heat capacity at first wall layer, 570 !< above ground floor level 571 INTEGER(iwp) :: ind_hc1_gfl = 25 !< index in input list for heat capacity at first wall layer, ground floor level 572 INTEGER(iwp) :: ind_hc1_wall_r = 89 !< index in input list for heat capacity at first wall layer, roof 573 INTEGER(iwp) :: ind_hc1_win_agfl = 71 !< index in input list for heat capacity at first window layer, 574 !< above ground floor level 575 INTEGER(iwp) :: ind_hc1_win_gfl = 38 !< index in input list for heat capacity at first window layer, 576 !< ground floor level 577 INTEGER(iwp) :: ind_hc1_win_r = 102 !< index in input list for heat capacity at first window layer, roof 578 INTEGER(iwp) :: ind_hc2_agfl = 59 !< index in input list for heat capacity at second wall layer, 579 !< above ground floor level 580 INTEGER(iwp) :: ind_hc2_gfl = 26 !< index in input list for heat capacity at second wall layer, ground floor level 581 INTEGER(iwp) :: ind_hc2_wall_r = 90 !< index in input list for heat capacity at second wall layer, roof 582 INTEGER(iwp) :: ind_hc2_win_agfl = 72 !< index in input list for heat capacity at second window layer, 583 !< above ground floor level 584 INTEGER(iwp) :: ind_hc2_win_gfl = 39 !< index in input list for heat capacity at second window layer, 585 !< ground floor level 586 INTEGER(iwp) :: ind_hc2_win_r = 103 !< index in input list for heat capacity at second window layer, roof 587 INTEGER(iwp) :: ind_hc3_agfl = 60 !< index in input list for heat capacity at third wall layer, 588 !< above ground floor level 589 INTEGER(iwp) :: ind_hc3_gfl = 27 !< index in input list for heat capacity at third wall layer, ground floor level 590 INTEGER(iwp) :: ind_hc3_wall_r = 91 !< index in input list for heat capacity at third wall layer, roof 591 INTEGER(iwp) :: ind_hc3_win_agfl = 73 !< index in input list for heat capacity at third window layer, 592 !< above ground floor level 593 INTEGER(iwp) :: ind_hc3_win_gfl = 40 !< index in input list for heat capacity at third window layer, 594 !< ground floor level 595 INTEGER(iwp) :: ind_hc3_win_r = 104 !< index in input list for heat capacity at third window layer, roof 596 INTEGER(iwp) :: ind_gflh = 17 !< index in input list for ground floor level height 597 INTEGER(iwp) :: ind_lai_r_agfl = 113 !< index in input list for LAI on roof, above ground floor level 584 598 INTEGER(iwp) :: ind_lai_r_gfl = 113 !< index in input list for LAI on roof, ground floor level 585 INTEGER(iwp) :: ind_lai_w_agfl = 81 !< index in input list for LAI on wall, above ground floor level 586 INTEGER(iwp) :: ind_lai_w_gfl = 48 !< index in input list for LAI on wall, ground floor level 587 INTEGER(iwp) :: ind_tc1_agfl = 61 !< index in input list for thermal conductivity at first wall layer, above ground floor level 588 INTEGER(iwp) :: ind_tc1_gfl = 28 !< index in input list for thermal conductivity at first wall layer, ground floor level 589 INTEGER(iwp) :: ind_tc1_wall_r = 92 !< index in input list for thermal conductivity at first wall layer, roof 590 INTEGER(iwp) :: ind_tc1_win_agfl = 74 !< index in input list for thermal conductivity at first window layer, above ground floor level 591 INTEGER(iwp) :: ind_tc1_win_gfl = 41 !< index in input list for thermal conductivity at first window layer, ground floor level 592 INTEGER(iwp) :: ind_tc1_win_r = 105 !< index in input list for thermal conductivity at first window layer, roof 593 INTEGER(iwp) :: ind_tc2_agfl = 62 !< index in input list for thermal conductivity at second wall layer, above ground floor level 594 INTEGER(iwp) :: ind_tc2_gfl = 29 !< index in input list for thermal conductivity at second wall layer, ground floor level 595 INTEGER(iwp) :: ind_tc2_wall_r = 93 !< index in input list for thermal conductivity at second wall layer, roof 596 INTEGER(iwp) :: ind_tc2_win_agfl = 75 !< index in input list for thermal conductivity at second window layer, above ground floor level 597 INTEGER(iwp) :: ind_tc2_win_gfl = 42 !< index in input list for thermal conductivity at second window layer, ground floor level 598 INTEGER(iwp) :: ind_tc2_win_r = 106 !< index in input list for thermal conductivity at second window layer, ground floor level 599 INTEGER(iwp) :: ind_tc3_agfl = 63 !< index in input list for thermal conductivity at third wall layer, above ground floor level 600 INTEGER(iwp) :: ind_tc3_gfl = 30 !< index in input list for thermal conductivity at third wall layer, ground floor level 601 INTEGER(iwp) :: ind_tc3_wall_r = 94 !< index in input list for thermal conductivity at third wall layer, roof 602 INTEGER(iwp) :: ind_tc3_win_agfl = 76 !< index in input list for thermal conductivity at third window layer, above ground floor level 603 INTEGER(iwp) :: ind_tc3_win_gfl = 43 !< index in input list for thermal conductivity at third window layer, ground floor level 604 INTEGER(iwp) :: ind_tc3_win_r = 107 !< index in input list for thermal conductivity at third window layer, roof 605 INTEGER(iwp) :: ind_thick_1_agfl = 54 !< index for wall layer thickness - 1st layer above ground floor level 606 INTEGER(iwp) :: ind_thick_1_gfl = 21 !< index for wall layer thickness - 1st layer ground floor level 607 INTEGER(iwp) :: ind_thick_1_wall_r = 85 !< index for wall layer thickness - 1st layer roof 608 INTEGER(iwp) :: ind_thick_1_win_agfl = 67 !< index for window layer thickness - 1st layer above ground floor level 609 INTEGER(iwp) :: ind_thick_1_win_gfl = 34 !< index for window layer thickness - 1st layer ground floor level 610 INTEGER(iwp) :: ind_thick_1_win_r = 98 !< index for window layer thickness - 1st layer roof 611 INTEGER(iwp) :: ind_thick_2_agfl = 55 !< index for wall layer thickness - 2nd layer above ground floor level 612 INTEGER(iwp) :: ind_thick_2_gfl = 22 !< index for wall layer thickness - 2nd layer ground floor level 613 INTEGER(iwp) :: ind_thick_2_wall_r = 86 !< index for wall layer thickness - 2nd layer roof 614 INTEGER(iwp) :: ind_thick_2_win_agfl = 68 !< index for window layer thickness - 2nd layer above ground floor level 615 INTEGER(iwp) :: ind_thick_2_win_gfl = 35 !< index for window layer thickness - 2nd layer ground floor level 616 INTEGER(iwp) :: ind_thick_2_win_r = 99 !< index for window layer thickness - 2nd layer roof 617 INTEGER(iwp) :: ind_thick_3_agfl = 56 !< index for wall layer thickness - 3rd layer above ground floor level 618 INTEGER(iwp) :: ind_thick_3_gfl = 23 !< index for wall layer thickness - 3rd layer ground floor level 619 INTEGER(iwp) :: ind_thick_3_wall_r = 87 !< index for wall layer thickness - 3rd layer roof 620 INTEGER(iwp) :: ind_thick_3_win_agfl = 69 !< index for window layer thickness - 3rd layer above ground floor level 621 INTEGER(iwp) :: ind_thick_3_win_gfl = 36 !< index for window layer thickness - 3rd layer ground floor level 622 INTEGER(iwp) :: ind_thick_3_win_r = 100 !< index for window layer thickness - 3rd layer roof 623 INTEGER(iwp) :: ind_thick_4_agfl = 57 !< index for wall layer thickness - 4th layer above ground floor level 624 INTEGER(iwp) :: ind_thick_4_gfl = 24 !< index for wall layer thickness - 4th layer ground floor level 625 INTEGER(iwp) :: ind_thick_4_wall_r = 88 !< index for wall layer thickness - 4st layer roof 626 INTEGER(iwp) :: ind_thick_4_win_agfl = 70 !< index for window layer thickness - 4th layer above ground floor level 627 INTEGER(iwp) :: ind_thick_4_win_gfl = 37 !< index for window layer thickness - 4th layer ground floor level 628 INTEGER(iwp) :: ind_thick_4_win_r = 101 !< index for window layer thickness - 4th layer roof 629 INTEGER(iwp) :: ind_trans_agfl = 78 !< index in input list for window transmissivity, above ground floor level 630 INTEGER(iwp) :: ind_trans_gfl = 45 !< index in input list for window transmissivity, ground floor level 631 INTEGER(iwp) :: ind_trans_r = 109 !< index in input list for window transmissivity, roof 632 INTEGER(iwp) :: ind_wall_frac_agfl = 53 !< index in input list for wall fraction, above ground floor level 633 INTEGER(iwp) :: ind_wall_frac_gfl = 20 !< index in input list for wall fraction, ground floor level 634 INTEGER(iwp) :: ind_wall_frac_r = 84 !< index in input list for wall fraction, roof 635 INTEGER(iwp) :: ind_win_frac_agfl = 66 !< index in input list for window fraction, above ground floor level 636 INTEGER(iwp) :: ind_win_frac_gfl = 33 !< index in input list for window fraction, ground floor level 637 INTEGER(iwp) :: ind_win_frac_r = 97 !< index in input list for window fraction, roof 638 INTEGER(iwp) :: ind_z0_agfl = 51 !< index in input list for z0, above ground floor level 639 INTEGER(iwp) :: ind_z0_gfl = 18 !< index in input list for z0, ground floor level 640 INTEGER(iwp) :: ind_z0qh_agfl = 52 !< index in input list for z0h / z0q, above ground floor level 641 INTEGER(iwp) :: ind_z0qh_gfl = 19 !< index in input list for z0h / z0q, ground floor level 642 INTEGER(iwp) :: ind_green_type_roof = 116 !< index in input list for type of green roof 643 644 645 REAL(wp) :: roof_height_limit = 4.0_wp !< height for distinguish between land surfaces and roofs 599 INTEGER(iwp) :: ind_lai_w_agfl = 81 !< index in input list for LAI on wall, above ground floor level 600 INTEGER(iwp) :: ind_lai_w_gfl = 48 !< index in input list for LAI on wall, ground floor level 601 INTEGER(iwp) :: ind_tc1_agfl = 61 !< index in input list for thermal conductivity at first wall layer, 602 !< above ground floor level 603 INTEGER(iwp) :: ind_tc1_gfl = 28 !< index in input list for thermal conductivity at first wall layer, 604 !< ground floor level 605 INTEGER(iwp) :: ind_tc1_wall_r = 92 !< index in input list for thermal conductivity at first wall layer, roof 606 INTEGER(iwp) :: ind_tc1_win_agfl = 74 !< index in input list for thermal conductivity at first window layer, 607 !< above ground floor level 608 INTEGER(iwp) :: ind_tc1_win_gfl = 41 !< index in input list for thermal conductivity at first window layer, 609 !< ground floor level 610 INTEGER(iwp) :: ind_tc1_win_r = 105 !< index in input list for thermal conductivity at first window layer, roof 611 INTEGER(iwp) :: ind_tc2_agfl = 62 !< index in input list for thermal conductivity at second wall layer, 612 !< above ground floor level 613 INTEGER(iwp) :: ind_tc2_gfl = 29 !< index in input list for thermal conductivity at second wall layer, 614 !< ground floor level 615 INTEGER(iwp) :: ind_tc2_wall_r = 93 !< index in input list for thermal conductivity at second wall layer, roof 616 INTEGER(iwp) :: ind_tc2_win_agfl = 75 !< index in input list for thermal conductivity at second window layer, 617 !< above ground floor level 618 INTEGER(iwp) :: ind_tc2_win_gfl = 42 !< index in input list for thermal conductivity at second window layer, 619 !< ground floor level 620 INTEGER(iwp) :: ind_tc2_win_r = 106 !< index in input list for thermal conductivity at second window layer, 621 !< ground floor level 622 INTEGER(iwp) :: ind_tc3_agfl = 63 !< index in input list for thermal conductivity at third wall layer, 623 !< above ground floor level 624 INTEGER(iwp) :: ind_tc3_gfl = 30 !< index in input list for thermal conductivity at third wall layer, 625 !< ground floor level 626 INTEGER(iwp) :: ind_tc3_wall_r = 94 !< index in input list for thermal conductivity at third wall layer, roof 627 INTEGER(iwp) :: ind_tc3_win_agfl = 76 !< index in input list for thermal conductivity at third window layer, 628 !< above ground floor level 629 INTEGER(iwp) :: ind_tc3_win_gfl = 43 !< index in input list for thermal conductivity at third window layer, 630 !< ground floor level 631 INTEGER(iwp) :: ind_tc3_win_r = 107 !< index in input list for thermal conductivity at third window layer, roof 632 INTEGER(iwp) :: ind_thick_1_agfl = 54 !< index for wall layer thickness - 1st layer above ground floor level 633 INTEGER(iwp) :: ind_thick_1_gfl = 21 !< index for wall layer thickness - 1st layer ground floor level 634 INTEGER(iwp) :: ind_thick_1_wall_r = 85 !< index for wall layer thickness - 1st layer roof 635 INTEGER(iwp) :: ind_thick_1_win_agfl = 67 !< index for window layer thickness - 1st layer above ground floor level 636 INTEGER(iwp) :: ind_thick_1_win_gfl = 34 !< index for window layer thickness - 1st layer ground floor level 637 INTEGER(iwp) :: ind_thick_1_win_r = 98 !< index for window layer thickness - 1st layer roof 638 INTEGER(iwp) :: ind_thick_2_agfl = 55 !< index for wall layer thickness - 2nd layer above ground floor level 639 INTEGER(iwp) :: ind_thick_2_gfl = 22 !< index for wall layer thickness - 2nd layer ground floor level 640 INTEGER(iwp) :: ind_thick_2_wall_r = 86 !< index for wall layer thickness - 2nd layer roof 641 INTEGER(iwp) :: ind_thick_2_win_agfl = 68 !< index for window layer thickness - 2nd layer above ground floor level 642 INTEGER(iwp) :: ind_thick_2_win_gfl = 35 !< index for window layer thickness - 2nd layer ground floor level 643 INTEGER(iwp) :: ind_thick_2_win_r = 99 !< index for window layer thickness - 2nd layer roof 644 INTEGER(iwp) :: ind_thick_3_agfl = 56 !< index for wall layer thickness - 3rd layer above ground floor level 645 INTEGER(iwp) :: ind_thick_3_gfl = 23 !< index for wall layer thickness - 3rd layer ground floor level 646 INTEGER(iwp) :: ind_thick_3_wall_r = 87 !< index for wall layer thickness - 3rd layer roof 647 INTEGER(iwp) :: ind_thick_3_win_agfl = 69 !< index for window layer thickness - 3rd layer above ground floor level 648 INTEGER(iwp) :: ind_thick_3_win_gfl = 36 !< index for window layer thickness - 3rd layer ground floor level 649 INTEGER(iwp) :: ind_thick_3_win_r = 100 !< index for window layer thickness - 3rd layer roof 650 INTEGER(iwp) :: ind_thick_4_agfl = 57 !< index for wall layer thickness - 4th layer above ground floor level 651 INTEGER(iwp) :: ind_thick_4_gfl = 24 !< index for wall layer thickness - 4th layer ground floor level 652 INTEGER(iwp) :: ind_thick_4_wall_r = 88 !< index for wall layer thickness - 4st layer roof 653 INTEGER(iwp) :: ind_thick_4_win_agfl = 70 !< index for window layer thickness - 4th layer above ground floor level 654 INTEGER(iwp) :: ind_thick_4_win_gfl = 37 !< index for window layer thickness - 4th layer ground floor level 655 INTEGER(iwp) :: ind_thick_4_win_r = 101 !< index for window layer thickness - 4th layer roof 656 INTEGER(iwp) :: ind_trans_agfl = 78 !< index in input list for window transmissivity, above ground floor level 657 INTEGER(iwp) :: ind_trans_gfl = 45 !< index in input list for window transmissivity, ground floor level 658 INTEGER(iwp) :: ind_trans_r = 109 !< index in input list for window transmissivity, roof 659 INTEGER(iwp) :: ind_wall_frac_agfl = 53 !< index in input list for wall fraction, above ground floor level 660 INTEGER(iwp) :: ind_wall_frac_gfl = 20 !< index in input list for wall fraction, ground floor level 661 INTEGER(iwp) :: ind_wall_frac_r = 84 !< index in input list for wall fraction, roof 662 INTEGER(iwp) :: ind_win_frac_agfl = 66 !< index in input list for window fraction, above ground floor level 663 INTEGER(iwp) :: ind_win_frac_gfl = 33 !< index in input list for window fraction, ground floor level 664 INTEGER(iwp) :: ind_win_frac_r = 97 !< index in input list for window fraction, roof 665 INTEGER(iwp) :: ind_z0_agfl = 51 !< index in input list for z0, above ground floor level 666 INTEGER(iwp) :: ind_z0_gfl = 18 !< index in input list for z0, ground floor level 667 INTEGER(iwp) :: ind_z0qh_agfl = 52 !< index in input list for z0h / z0q, above ground floor level 668 INTEGER(iwp) :: ind_z0qh_gfl = 19 !< index in input list for z0h / z0q, ground floor level 669 INTEGER(iwp) :: ind_green_type_roof = 116 !< index in input list for type of green roof 670 671 672 REAL(wp) :: roof_height_limit = 4.0_wp !< height for distinguish between land surfaces and roofs 646 673 REAL(wp) :: ground_floor_level = 4.0_wp !< default ground floor level 647 674 648 675 649 676 CHARACTER(37), DIMENSION(0:7), PARAMETER :: building_type_name = (/ & 650 'user-defined ', & !0651 'residential - 1950 ', & !1652 'residential 1951 - 2000 ', & !2653 'residential 2001 - ', & !3654 'office - 1950 ', & !4655 'office 1951 - 2000 ', & !5656 'office 2001 - ', & !6657 'bridges ' & !7677 'user-defined ', & !< type 0 678 'residential - 1950 ', & !< type 1 679 'residential 1951 - 2000 ', & !< type 2 680 'residential 2001 - ', & !< type 3 681 'office - 1950 ', & !< type 4 682 'office 1951 - 2000 ', & !< type 5 683 'office 2001 - ', & !< type 6 684 'bridges ' & !< type 7 658 685 /) 659 686 ! … … 744 771 0.57_wp, 0.57_wp, 0.57_wp, 0.91_wp, & !parameter 41-44 745 772 0.75_wp, 27.0_wp, 0.0_wp, 1.5_wp, 0.86_wp, & !parameter 45-49 746 5.0_wp, 0.001_wp, 0.0001_wp, 0.7_wp, 0.005_wp, & !parameter 50-54773 5.0_wp, 0.001_wp, 0.0001_wp, 0.7_wp, 0.005_wp, & !parameter 50-54 747 774 0.01_wp, 0.39_wp, 0.63_wp, 2200000.0_wp, & !parameter 55-58 748 775 1400000.0_wp, 1300000.0_wp, 0.35_wp, 0.8_wp, & !parameter 59-62 … … 757 784 1736000.0_wp, 1736000.0_wp, 0.57_wp, 0.57_wp, 0.57_wp, & !parameter 103-107 758 785 0.91_wp, 0.75_wp, 27.0_wp, 0.0_wp, 0.0_wp, 1.5_wp, & !parameter 108-113 759 0.86_wp, 5.0_wp, 0.0_wp, & !parameter 114-116786 0.86_wp, 5.0_wp, 0.0_wp, & !parameter 114-116 760 787 299.15_wp, 293.15_wp, 0.8_wp, 0.76_wp, 5.0_wp, & !parameter 117-121 761 788 0.1_wp, 0.5_wp, 0.0_wp, 3.5_wp, 370000.0_wp, 4.5_wp, & !parameter 122-127 … … 772 799 0.11_wp, 0.11_wp, 0.11_wp, 0.11_wp, & !parameter 41-44 773 800 0.7_wp, 27.0_wp, 0.0_wp, 1.5_wp, 0.86_wp, & !parameter 45-49 774 5.0_wp, 0.001_wp, 0.0001_wp, 0.73_wp, 0.005_wp, & !parameter 50-54801 5.0_wp, 0.001_wp, 0.0001_wp, 0.73_wp, 0.005_wp, & !parameter 50-54 775 802 0.01_wp, 0.31_wp, 0.43_wp, 2000000.0_wp, & !parameter 55-58 776 803 103000.0_wp, 900000.0_wp, 0.35_wp, 0.38_wp, & !parameter 59-62 … … 785 812 1736000.0_wp, 1736000.0_wp, 0.11_wp, 0.11_wp, 0.11_wp, & !parameter 103-107 786 813 0.87_wp, 0.7_wp, 27.0_wp, 0.0_wp, 0.0_wp, 1.5_wp, & !parameter 108-113 787 0.86_wp, 5.0_wp, 0.0_wp, & !parameter 114-116814 0.86_wp, 5.0_wp, 0.0_wp, & !parameter 114-116 788 815 299.15_wp, 293.15_wp, 0.8_wp, 0.6_wp, 3.0_wp, & !parameter 117-121 789 816 0.1_wp, 0.5_wp, 0.0_wp, 2.5_wp, 165000.0_wp, 4.5_wp, & !parameter 122-127 … … 798 825 27.0_wp, 0.25_wp, 0.003_wp, 0.006_wp, 0.012_wp, 0.018_wp, & !parameter 32-37 799 826 1736000.0_wp, 1736000.0_wp, 1736000.0_wp, & !parameter 38-40 800 0.037_wp, 0.037_wp, 0.037_wp, 0.8_wp, 827 0.037_wp, 0.037_wp, 0.037_wp, 0.8_wp, & !parameter 41-44 801 828 0.6_wp, 27.0_wp, 0.0_wp, 1.5_wp, 0.86_wp, & !parameter 45-49 802 5.0_wp, 0.001_wp, 0.0001_wp, 0.7_wp, 0.005_wp, & !parameter 50-54829 5.0_wp, 0.001_wp, 0.0001_wp, 0.7_wp, 0.005_wp, & !parameter 50-54 803 830 0.01_wp, 0.41_wp, 0.7_wp, 2000000.0_wp, & !parameter 55-58 804 831 103000.0_wp, 900000.0_wp, 0.35_wp, 0.14_wp, & !parameter 59-62 … … 806 833 0.003_wp, 0.006_wp, 0.012_wp, 0.018_wp, & !parameter 67-70 807 834 1736000.0_wp, 1736000.0_wp, 1736000.0_wp, & !parameter 71-73 808 0.037_wp, 0.037_wp, 0.037_wp, 0.8_wp, 0.6_wp, 835 0.037_wp, 0.037_wp, 0.037_wp, 0.8_wp, 0.6_wp, & !parameter 74-78 809 836 27.0_wp, 0.0_wp, 1.5_wp, 0.86_wp, 5.0_wp, 1.0_wp, & !parameter 79-84 810 837 0.005_wp, 0.01_wp, 0.41_wp, 0.7_wp, 2000000.0_wp, 103000.0_wp, & !parameter 85-90 811 838 900000.0_wp, 0.35_wp, 0.14_wp, 0.035_wp, 0.93_wp, 27.0_wp, 0.0_wp, & !parameter 91-97 812 839 0.003_wp, 0.006_wp, 0.012_wp, 0.018_wp, 1736000.0_wp, & !parameter 98-102 813 1736000.0_wp, 1736000.0_wp, 0.037_wp, 0.037_wp, 0.037_wp, 840 1736000.0_wp, 1736000.0_wp, 0.037_wp, 0.037_wp, 0.037_wp, & !parameter 103-107 814 841 0.8_wp, 0.6_wp, 27.0_wp, 0.0_wp, 0.0_wp, 1.5_wp, & !parameter 108-113 815 0.86_wp, 5.0_wp, 0.0_wp, & !parameter 114-116842 0.86_wp, 5.0_wp, 0.0_wp, & !parameter 114-116 816 843 299.15_wp, 293.15_wp, 0.8_wp, 0.5_wp, 0.6_wp, & !parameter 117-121 817 844 0.1_wp, 0.5_wp, 0.8_wp, 2.5_wp, 80000.0_wp, 4.5_wp, & !parameter 122-127 … … 828 855 0.57_wp, 0.57_wp, 0.57_wp, 0.91_wp, & !parameter 41-44 829 856 0.75_wp, 27.0_wp, 0.0_wp, 1.5_wp, 0.86_wp, & !parameter 45-49 830 5.0_wp, 0.001_wp, 0.0001_wp, 0.5_wp, 0.005_wp, & !parameter 50-54857 5.0_wp, 0.001_wp, 0.0001_wp, 0.5_wp, 0.005_wp, & !parameter 50-54 831 858 0.01_wp, 0.39_wp, 0.63_wp, 2200000.0_wp, & !parameter 55-58 832 859 1400000.0_wp, 1300000.0_wp, 0.35_wp, 0.8_wp, & !parameter 59-62 … … 841 868 1736000.0_wp, 1736000.0_wp, 0.57_wp, 0.57_wp, 0.57_wp, & !parameter 103-107 842 869 0.91_wp, 0.75_wp, 27.0_wp, 0.0_wp, 0.0_wp, 1.5_wp, & !parameter 108-113 843 0.86_wp, 5.0_wp, 0.0_wp, & !parameter 114-116870 0.86_wp, 5.0_wp, 0.0_wp, & !parameter 114-116 844 871 299.15_wp, 293.15_wp, 0.8_wp, 0.76_wp, 5.0_wp, & !parameter 117-121 845 872 0.1_wp, 1.5_wp, 0.0_wp, 3.5_wp, 370000.0_wp, 4.5_wp, & !parameter 122-127 … … 856 883 0.11_wp, 0.11_wp, 0.11_wp, 0.87_wp, & !parameter 41-44 857 884 0.7_wp, 27.0_wp, 0.0_wp, 1.5_wp, 0.86_wp, & !parameter 45-49 858 5.0_wp, 0.001_wp, 0.0001_wp, 0.5_wp, 0.005_wp, & !parameter 50-54885 5.0_wp, 0.001_wp, 0.0001_wp, 0.5_wp, 0.005_wp, & !parameter 50-54 859 886 0.01_wp, 0.31_wp, 0.43_wp, 2000000.0_wp, & !parameter 55-58 860 887 103000.0_wp, 900000.0_wp, 0.35_wp, 0.38_wp, & !parameter 59-62 … … 869 896 1736000.0_wp, 1736000.0_wp, 0.11_wp, 0.11_wp, 0.11_wp, & !parameter 103-107 870 897 0.87_wp, 0.7_wp, 27.0_wp, 0.0_wp, 0.0_wp, 1.5_wp, & !parameter 108-113 871 0.86_wp, 5.0_wp, 0.0_wp, & !parameter 114-116898 0.86_wp, 5.0_wp, 0.0_wp, & !parameter 114-116 872 899 299.15_wp, 293.15_wp, 0.8_wp, 0.6_wp, 3.0_wp, & !parameter 117-121 873 900 0.1_wp, 1.5_wp, 0.65_wp, 2.5_wp, 165000.0_wp, 4.5_wp, & !parameter 122-127 … … 882 909 27.0_wp, 0.525_wp, 0.003_wp, 0.006_wp, 0.012_wp, 0.018_wp, & !parameter 32-37 883 910 1736000.0_wp, 1736000.0_wp, 1736000.0_wp, & !parameter 38-40 884 0.037_wp, 0.037_wp, 0.037_wp, 0.8_wp, 911 0.037_wp, 0.037_wp, 0.037_wp, 0.8_wp, & !parameter 41-44 885 912 0.6_wp, 27.0_wp, 0.0_wp, 1.5_wp, 0.86_wp, & !parameter 45-49 886 5.0_wp, 0.001_wp, 0.0001_wp, 0.425_wp, 0.005_wp, & !parameter 50-54913 5.0_wp, 0.001_wp, 0.0001_wp, 0.425_wp, 0.005_wp, & !parameter 50-54 887 914 0.01_wp, 0.41_wp, 0.7_wp, 2000000.0_wp, & !parameter 55-58 888 915 103000.0_wp, 900000.0_wp, 0.35_wp, 0.14_wp, & !parameter 59-62 … … 890 917 0.003_wp, 0.006_wp, 0.012_wp, 0.018_wp, & !parameter 67-70 891 918 1736000.0_wp, 1736000.0_wp, 1736000.0_wp, & !parameter 71-73 892 0.037_wp, 0.037_wp, 0.037_wp, 0.8_wp, 0.6_wp, 919 0.037_wp, 0.037_wp, 0.037_wp, 0.8_wp, 0.6_wp, & !parameter 74-78 893 920 27.0_wp, 0.0_wp, 1.5_wp, 0.86_wp, 5.0_wp, 1.0_wp, & !parameter 79-84 894 921 0.005_wp, 0.01_wp, 0.41_wp, 0.7_wp, 2000000.0_wp, 103000.0_wp, & !parameter 85-90 895 922 900000.0_wp, 0.35_wp, 0.14_wp, 0.035_wp, 0.91_wp, 27.0_wp, 0.0_wp, & !parameter 91-97 896 923 0.003_wp, 0.006_wp, 0.012_wp, 0.018_wp, 1736000.0_wp, & !parameter 98-102 897 1736000.0_wp, 1736000.0_wp, 0.037_wp, 0.037_wp, 0.037_wp, 924 1736000.0_wp, 1736000.0_wp, 0.037_wp, 0.037_wp, 0.037_wp, & !parameter 103-107 898 925 0.8_wp, 0.6_wp, 27.0_wp, 0.0_wp, 0.0_wp, 1.5_wp, & !parameter 108-113 899 0.86_wp, 5.0_wp, 0.0_wp, & !parameter 114-116926 0.86_wp, 5.0_wp, 0.0_wp, & !parameter 114-116 900 927 299.15_wp, 293.15_wp, 0.8_wp, 0.5_wp, 0.6_wp, & !parameter 117-121 901 928 0.1_wp, 1.5_wp, 0.9_wp, 2.5_wp, 80000.0_wp, 4.5_wp, & !parameter 122-127 … … 912 939 0.57_wp, 0.57_wp, 0.57_wp, 0.8_wp, & !parameter 41-44 913 940 0.6_wp, 27.0_wp, 0.0_wp, 1.5_wp, 0.86_wp, & !parameter 45-49 914 5.0_wp, 0.001_wp, 0.0001_wp, 1.0_wp, 0.29_wp, & !parameter 50-54941 5.0_wp, 0.001_wp, 0.0001_wp, 1.0_wp, 0.29_wp, & !parameter 50-54 915 942 0.295_wp, 0.695_wp, 0.985_wp, 1950400.0_wp, & !parameter 55-58 916 943 1848000.0_wp, 1848000.0_wp, 0.7_wp, 1.0_wp, & !parameter 59-62 … … 925 952 1736000.0_wp, 1736000.0_wp, 0.57_wp, 0.57_wp, 0.57_wp, & !parameter 103-107 926 953 0.8_wp, 0.6_wp, 27.0_wp, 0.0_wp, 0.0_wp, 1.5_wp, & !parameter 108-113 927 0.86_wp, 5.0_wp, 0.0_wp, & !parameter 114-116954 0.86_wp, 5.0_wp, 0.0_wp, & !parameter 114-116 928 955 299.15_wp, 293.15_wp, 0.8_wp, 100.0_wp, 100.0_wp, & !parameter 117-121 929 956 20.0_wp, 20.0_wp, 0.0_wp, 1.0_wp, 1.0_wp, 4.5_wp, & !parameter 122-127 … … 944 971 945 972 TYPE surf_type_usm 946 REAL(wp), DIMENSION(:), ALLOCATABLE :: var_usm_1d !< 1D prognostic variable947 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: var_usm_2d !< 2D prognostic variable973 REAL(wp), DIMENSION(:), ALLOCATABLE :: var_usm_1d !< 1D prognostic variable 974 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: var_usm_2d !< 2D prognostic variable 948 975 END TYPE surf_type_usm 949 976 950 TYPE(surf_type_usm), POINTER :: m_liq_usm_h, & !< liquid water reservoir (m), horizontal surface elements951 m_liq_usm_h_p !< progn. liquid water reservoir (m), horizontal surface elements952 953 TYPE(surf_type_usm), TARGET :: m_liq_usm_h_1, & !<954 m_liq_usm_h_2 !<955 956 TYPE(surf_type_usm), DIMENSION(:), POINTER :: &957 m_liq_usm_v, & !< liquid water reservoir (m), vertical surface elements958 m_liq_usm_v_p !< progn. liquid water reservoir (m), vertical surface elements959 960 TYPE(surf_type_usm), DIMENSION(0:3), TARGET :: &961 m_liq_usm_v_1, & !<962 m_liq_usm_v_2 !<977 TYPE(surf_type_usm), POINTER :: m_liq_usm_h, & !< liquid water reservoir (m), horizontal surface elements 978 m_liq_usm_h_p !< progn. liquid water reservoir (m), horizontal surface elements 979 980 TYPE(surf_type_usm), TARGET :: m_liq_usm_h_1, & !< 981 m_liq_usm_h_2 !< 982 983 TYPE(surf_type_usm), DIMENSION(:), POINTER :: & 984 m_liq_usm_v, & !< liquid water reservoir (m), vertical surface elements 985 m_liq_usm_v_p !< progn. liquid water reservoir (m), vertical surface elements 986 987 TYPE(surf_type_usm), DIMENSION(0:3), TARGET :: & 988 m_liq_usm_v_1, & !< 989 m_liq_usm_v_2 !< 963 990 964 991 TYPE(surf_type_usm), TARGET :: tm_liq_usm_h_m !< liquid water reservoir tendency (m), horizontal surface elements 965 TYPE(surf_type_usm), DIMENSION(0:3), TARGET :: tm_liq_usm_v_m !< liquid water reservoir tendency (m), vertical surface elements966 967 968 ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!992 TYPE(surf_type_usm), DIMENSION(0:3), TARGET :: tm_liq_usm_v_m !< liquid water reservoir tendency (m), 993 !< vertical surface elements 994 995 ! 969 996 !-- anthropogenic heat sources 970 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!971 997 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: aheat !< daily average of anthropogenic heat (W/m2) 972 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: aheatprof !< diurnal profiles of anthropogenic heat for particular layers 998 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: aheatprof !< diurnal profiles of anthropogenic heat 999 !< for particular layers 973 1000 INTEGER(iwp) :: naheatlayers = 1 !< number of layers of anthropogenic heat 974 1001 975 ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!1002 ! 976 1003 !-- wall surface model 977 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!978 1004 !-- wall surface model constants 979 1005 INTEGER(iwp), PARAMETER :: nzb_wall = 0 !< inner side of the wall model (to be switched) … … 981 1007 INTEGER(iwp), PARAMETER :: nzw = 4 !< number of wall layers (fixed for now) 982 1008 983 REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default = (/0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp /) 984 !< normalized soil, wall and roof layer depths (m/m) 985 ! REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default = (/0.33_wp, 0.66_wp, 1.0_wp /) 986 REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default_window = (/0.25_wp, 0.5_wp, 0.75_wp, 1.0_wp /) 987 ! REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default_window = (/0.33_wp, 0.66_wp, 1.0_wp /) 988 ! REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default_window = (/0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp /) 989 !< normalized window layer depths (m/m) 990 ! REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default_green = (/0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp /) 991 !< normalized green layer depths (m/m) 992 REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default_green = (/0.25_wp, 0.5_wp, 0.75_wp, 1.0_wp /) 993 ! REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default_green = (/0.33_wp, 0.66_wp, 1.0_wp /) 994 995 996 REAL(wp) :: wall_inner_temperature = 295.0_wp !< temperature of the inner wall surface (~22 degrees C) (K) 997 REAL(wp) :: roof_inner_temperature = 295.0_wp !< temperature of the inner roof surface (~22 degrees C) (K) 998 REAL(wp) :: soil_inner_temperature = 288.0_wp !< temperature of the deep soil (~15 degrees C) (K) 999 REAL(wp) :: window_inner_temperature = 295.0_wp !< temperature of the inner window surface (~22 degrees C) (K) 1000 1001 REAL(wp) :: m_total = 0.0_wp !< weighted total water content of the soil (m3/m3) 1002 INTEGER(iwp) :: soil_type 1003 1004 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1009 REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default = (/0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp /) 1010 REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default_window = (/0.25_wp, 0.5_wp, 0.75_wp, 1.0_wp /) 1011 REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default_green = (/0.25_wp, 0.5_wp, 0.75_wp, 1.0_wp /) 1012 !< normalized soil, wall and roof, window and 1013 !<green layer depths (m/m) 1014 1015 REAL(wp) :: wall_inner_temperature = 295.0_wp !< temperature of the inner wall 1016 !< surface (~22 degrees C) (K) 1017 REAL(wp) :: roof_inner_temperature = 295.0_wp !< temperature of the inner roof 1018 !< surface (~22 degrees C) (K) 1019 REAL(wp) :: soil_inner_temperature = 288.0_wp !< temperature of the deep soil 1020 !< (~15 degrees C) (K) 1021 REAL(wp) :: window_inner_temperature = 295.0_wp !< temperature of the inner window 1022 !< surface (~22 degrees C) (K) 1023 1024 REAL(wp) :: m_total = 0.0_wp !< weighted total water content of the soil (m3/m3) 1025 INTEGER(iwp) :: soil_type 1026 1027 ! 1005 1028 !-- surface and material model variables for walls, ground, roofs 1006 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!1007 1029 REAL(wp), DIMENSION(:), ALLOCATABLE :: zwn !< normalized wall layer depths (m) 1008 1030 REAL(wp), DIMENSION(:), ALLOCATABLE :: zwn_window !< normalized window layer depths (m) … … 1023 1045 REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_green_h_2 1024 1046 1025 TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_wall_v1026 TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_wall_v_p1027 TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_window_v1028 TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_window_v_p1029 TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_green_v1030 TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_green_v_p1047 TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_wall_v 1048 TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_wall_v_p 1049 TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_window_v 1050 TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_window_v_p 1051 TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_green_v 1052 TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_green_v_p 1031 1053 1032 1054 TYPE(t_surf_vertical), DIMENSION(0:3), TARGET :: t_surf_wall_v_1 … … 1037 1059 TYPE(t_surf_vertical), DIMENSION(0:3), TARGET :: t_surf_green_v_2 1038 1060 1039 ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!1061 ! 1040 1062 !-- Energy balance variables 1041 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!1042 1063 !-- parameters of the land, roof and wall surfaces 1043 1064 … … 1062 1083 TYPE(t_wall_vertical), DIMENSION(0:3), TARGET :: swc_v_1, swc_v_2 1063 1084 1085 ! 1064 1086 !-- Surface and material parameters classes (surface_type) 1065 1087 !-- albedo, emissivity, lambda_surf, roughness, thickness, volumetric heat capacity, thermal conductivity 1066 INTEGER(iwp) :: n_surface_types !< number of the wall type categories 1067 INTEGER(iwp), PARAMETER :: n_surface_params = 9 !< number of parameters for each type of the wall 1068 INTEGER(iwp), PARAMETER :: ialbedo = 1 !< albedo of the surface 1069 INTEGER(iwp), PARAMETER :: iemiss = 2 !< emissivity of the surface 1070 INTEGER(iwp), PARAMETER :: ilambdas = 3 !< heat conductivity lambda S between surface and material ( W m-2 K-1 ) 1071 INTEGER(iwp), PARAMETER :: irough = 4 !< roughness length z0 for movements 1072 INTEGER(iwp), PARAMETER :: iroughh = 5 !< roughness length z0h for scalars (heat, humidity,...) 1073 INTEGER(iwp), PARAMETER :: icsurf = 6 !< Surface skin layer heat capacity (J m-2 K-1 ) 1074 INTEGER(iwp), PARAMETER :: ithick = 7 !< thickness of the surface (wall, roof, land) ( m ) 1075 INTEGER(iwp), PARAMETER :: irhoC = 8 !< volumetric heat capacity rho*C of the material ( J m-3 K-1 ) 1076 INTEGER(iwp), PARAMETER :: ilambdah = 9 !< thermal conductivity lambda H of the wall (W m-1 K-1 ) 1077 CHARACTER(12), DIMENSION(:), ALLOCATABLE :: surface_type_names !< names of wall types (used only for reports) 1078 INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: surface_type_codes !< codes of wall types 1079 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: surface_params !< parameters of wall types 1080 1081 1088 INTEGER(iwp) :: n_surface_types !< number of the wall type categories 1089 INTEGER(iwp), PARAMETER :: n_surface_params = 9 !< number of parameters for each type of the wall 1090 INTEGER(iwp), PARAMETER :: ialbedo = 1 !< albedo of the surface 1091 INTEGER(iwp), PARAMETER :: iemiss = 2 !< emissivity of the surface 1092 INTEGER(iwp), PARAMETER :: ilambdas = 3 !< heat conductivity lambda S between surface 1093 !< and material ( W m-2 K-1 ) 1094 INTEGER(iwp), PARAMETER :: irough = 4 !< roughness length z0 for movements 1095 INTEGER(iwp), PARAMETER :: iroughh = 5 !< roughness length z0h for scalars 1096 !< (heat, humidity,...) 1097 INTEGER(iwp), PARAMETER :: icsurf = 6 !< Surface skin layer heat capacity (J m-2 K-1 ) 1098 INTEGER(iwp), PARAMETER :: ithick = 7 !< thickness of the surface (wall, roof, land) ( m ) 1099 INTEGER(iwp), PARAMETER :: irhoC = 8 !< volumetric heat capacity rho*C of 1100 !< the material ( J m-3 K-1 ) 1101 INTEGER(iwp), PARAMETER :: ilambdah = 9 !< thermal conductivity lambda H 1102 !< of the wall (W m-1 K-1 ) 1103 CHARACTER(12), DIMENSION(:), ALLOCATABLE :: surface_type_names !< names of wall types (used only for reports) 1104 INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: surface_type_codes !< codes of wall types 1105 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: surface_params !< parameters of wall types 1106 1107 ! 1082 1108 !-- interfaces of subroutines accessed from outside of this module 1109 INTERFACE usm_3d_data_averaging 1110 MODULE PROCEDURE usm_3d_data_averaging 1111 END INTERFACE usm_3d_data_averaging 1112 1083 1113 INTERFACE usm_boundary_condition 1084 1114 MODULE PROCEDURE usm_boundary_condition … … 1105 1135 END INTERFACE usm_init 1106 1136 1137 INTERFACE usm_init_arrays 1138 MODULE PROCEDURE usm_init_arrays 1139 END INTERFACE usm_init_arrays 1140 1107 1141 INTERFACE usm_material_heat_model 1108 1142 MODULE PROCEDURE usm_material_heat_model … … 1133 1167 END INTERFACE usm_wrd_local 1134 1168 1135 INTERFACE usm_init_arrays1136 MODULE PROCEDURE usm_init_arrays1137 END INTERFACE usm_init_arrays1138 1139 INTERFACE usm_3d_data_averaging1140 MODULE PROCEDURE usm_3d_data_averaging1141 END INTERFACE usm_3d_data_averaging1142 1143 1169 1144 1170 SAVE 1145 1171 1146 1172 PRIVATE 1147 1173 1174 ! 1148 1175 !-- Public functions 1149 PUBLIC usm_boundary_condition, usm_check_parameters, usm_init, &1150 usm_rrd_local, &1151 usm_surface_energy_balance, usm_material_heat_model, &1152 usm_swap_timelevel, usm_check_data_output, usm_3d_data_averaging, &1153 usm_data_output_3d, usm_define_netcdf_grid, usm_parin, &1176 PUBLIC usm_boundary_condition, usm_check_parameters, usm_init, & 1177 usm_rrd_local, & 1178 usm_surface_energy_balance, usm_material_heat_model, & 1179 usm_swap_timelevel, usm_check_data_output, usm_3d_data_averaging, & 1180 usm_data_output_3d, usm_define_netcdf_grid, usm_parin, & 1154 1181 usm_wrd_local, usm_init_arrays 1155 1182 1183 ! 1156 1184 !-- Public parameters, constants and initial values 1157 PUBLIC usm_anthropogenic_heat, usm_material_model, usm_wall_mod, 1158 usm_green_heat_model, building_pars, 1159 nzb_wall, nzt_wall, t_wall_h, t_wall_v, 1185 PUBLIC usm_anthropogenic_heat, usm_material_model, usm_wall_mod, & 1186 usm_green_heat_model, building_pars, & 1187 nzb_wall, nzt_wall, t_wall_h, t_wall_v, & 1160 1188 t_window_h, t_window_v, building_type 1161 1189 … … 1181 1209 !-- Allocate radiation arrays which are part of the new data type. 1182 1210 !-- For horizontal surfaces. 1183 ALLOCATE ( surf_usm_h%surfhf(1:surf_usm_h%ns) )1184 ALLOCATE ( surf_usm_h%rad_net_l(1:surf_usm_h%ns) )1211 ALLOCATE ( surf_usm_h%surfhf(1:surf_usm_h%ns) ) 1212 ALLOCATE ( surf_usm_h%rad_net_l(1:surf_usm_h%ns) ) 1185 1213 ! 1186 1214 !-- For vertical surfaces 1187 1215 DO l = 0, 3 1188 ALLOCATE ( surf_usm_v(l)%surfhf(1:surf_usm_v(l)%ns) )1189 ALLOCATE ( surf_usm_v(l)%rad_net_l(1:surf_usm_v(l)%ns) )1216 ALLOCATE ( surf_usm_v(l)%surfhf(1:surf_usm_v(l)%ns) ) 1217 ALLOCATE ( surf_usm_v(l)%rad_net_l(1:surf_usm_v(l)%ns) ) 1190 1218 ENDDO 1191 1219 1220 ! 1192 1221 !-- Wall surface model 1193 1222 !-- allocate arrays for wall surface model and define pointers 1194 1195 1223 !-- allocate array of wall types and wall parameters 1196 1224 ALLOCATE ( surf_usm_h%surface_types(1:surf_usm_h%ns) ) … … 1200 1228 surf_usm_h%building_type_name = 'none' 1201 1229 DO l = 0, 3 1202 ALLOCATE ( surf_usm_v(l)%surface_types(1:surf_usm_v(l)%ns))1230 ALLOCATE ( surf_usm_v(l)%surface_types(1:surf_usm_v(l)%ns) ) 1203 1231 ALLOCATE ( surf_usm_v(l)%building_type(1:surf_usm_v(l)%ns) ) 1204 1232 ALLOCATE ( surf_usm_v(l)%building_type_name(1:surf_usm_v(l)%ns) ) … … 1209 1237 !-- Allocate albedo_type and albedo. Each surface element 1210 1238 !-- has 3 values, 0: wall fraction, 1: green fraction, 2: window fraction. 1211 ALLOCATE ( surf_usm_h%albedo_type(0:2,1:surf_usm_h%ns) )1212 ALLOCATE ( surf_usm_h%albedo(0:2,1:surf_usm_h%ns) )1239 ALLOCATE ( surf_usm_h%albedo_type(0:2,1:surf_usm_h%ns) ) 1240 ALLOCATE ( surf_usm_h%albedo(0:2,1:surf_usm_h%ns) ) 1213 1241 surf_usm_h%albedo_type = albedo_type 1214 1242 DO l = 0, 3 1215 ALLOCATE ( surf_usm_v(l)%albedo_type(0:2,1:surf_usm_v(l)%ns) )1216 ALLOCATE ( surf_usm_v(l)%albedo(0:2,1:surf_usm_v(l)%ns) )1243 ALLOCATE ( surf_usm_v(l)%albedo_type(0:2,1:surf_usm_v(l)%ns) ) 1244 ALLOCATE ( surf_usm_v(l)%albedo(0:2,1:surf_usm_v(l)%ns) ) 1217 1245 surf_usm_v(l)%albedo_type = albedo_type 1218 1246 ENDDO 1219 1247 1220 1221 1248 ! 1222 1249 !-- Allocate indoor target temperature for summer and winter 1223 ALLOCATE ( surf_usm_h%target_temp_summer(1:surf_usm_h%ns) )1224 ALLOCATE ( surf_usm_h%target_temp_winter(1:surf_usm_h%ns) )1250 ALLOCATE ( surf_usm_h%target_temp_summer(1:surf_usm_h%ns) ) 1251 ALLOCATE ( surf_usm_h%target_temp_winter(1:surf_usm_h%ns) ) 1225 1252 DO l = 0, 3 1226 ALLOCATE ( surf_usm_v(l)%target_temp_summer(1:surf_usm_v(l)%ns) )1227 ALLOCATE ( surf_usm_v(l)%target_temp_winter(1:surf_usm_v(l)%ns) )1253 ALLOCATE ( surf_usm_v(l)%target_temp_summer(1:surf_usm_v(l)%ns) ) 1254 ALLOCATE ( surf_usm_v(l)%target_temp_winter(1:surf_usm_v(l)%ns) ) 1228 1255 ENDDO 1229 1256 ! … … 1231 1258 ALLOCATE ( surf_usm_h%ground_level(1:surf_usm_h%ns) ) 1232 1259 DO l = 0, 3 1233 ALLOCATE ( surf_usm_v(l)%ground_level(1:surf_usm_v(l)%ns) )1260 ALLOCATE ( surf_usm_v(l)%ground_level(1:surf_usm_v(l)%ns) ) 1234 1261 ENDDO 1235 1262 ! 1236 1263 !-- Allocate arrays for relative surface fraction. 1237 1264 !-- 0 - wall fraction, 1 - green fraction, 2 - window fraction 1238 ALLOCATE ( surf_usm_h%frac(0:2,1:surf_usm_h%ns) )1265 ALLOCATE ( surf_usm_h%frac(0:2,1:surf_usm_h%ns) ) 1239 1266 surf_usm_h%frac = 0.0_wp 1240 1267 DO l = 0, 3 1241 ALLOCATE ( surf_usm_v(l)%frac(0:2,1:surf_usm_v(l)%ns) )1268 ALLOCATE ( surf_usm_v(l)%frac(0:2,1:surf_usm_v(l)%ns) ) 1242 1269 surf_usm_v(l)%frac = 0.0_wp 1243 1270 ENDDO 1244 1271 1272 ! 1245 1273 !-- wall and roof surface parameters. First for horizontal surfaces 1246 1274 ALLOCATE ( surf_usm_h%isroof_surf(1:surf_usm_h%ns) ) … … 1263 1291 !-- For vertical surfaces. 1264 1292 DO l = 0, 3 1265 ALLOCATE ( surf_usm_v(l)%lambda_surf(1:surf_usm_v(l)%ns) )1266 ALLOCATE ( surf_usm_v(l)%c_surface(1:surf_usm_v(l)%ns) )1293 ALLOCATE ( surf_usm_v(l)%lambda_surf(1:surf_usm_v(l)%ns) ) 1294 ALLOCATE ( surf_usm_v(l)%c_surface(1:surf_usm_v(l)%ns) ) 1267 1295 ALLOCATE ( surf_usm_v(l)%lambda_surf_window(1:surf_usm_v(l)%ns) ) 1268 1296 ALLOCATE ( surf_usm_v(l)%c_surface_window(1:surf_usm_v(l)%ns) ) 1269 1297 ALLOCATE ( surf_usm_v(l)%lambda_surf_green(1:surf_usm_v(l)%ns) ) 1270 1298 ALLOCATE ( surf_usm_v(l)%c_surface_green(1:surf_usm_v(l)%ns) ) 1271 ALLOCATE ( surf_usm_v(l)%transmissivity(1:surf_usm_v(l)%ns) )1272 ALLOCATE ( surf_usm_v(l)%lai(1:surf_usm_v(l)%ns) )1273 ALLOCATE ( surf_usm_v(l)%emissivity(0:2,1:surf_usm_v(l)%ns) )1274 ALLOCATE ( surf_usm_v(l)%r_a(1:surf_usm_v(l)%ns) )1275 ALLOCATE ( surf_usm_v(l)%r_a_green(1:surf_usm_v(l)%ns) )1276 ALLOCATE ( surf_usm_v(l)%r_a_window(1:surf_usm_v(l)%ns) )1299 ALLOCATE ( surf_usm_v(l)%transmissivity(1:surf_usm_v(l)%ns) ) 1300 ALLOCATE ( surf_usm_v(l)%lai(1:surf_usm_v(l)%ns) ) 1301 ALLOCATE ( surf_usm_v(l)%emissivity(0:2,1:surf_usm_v(l)%ns) ) 1302 ALLOCATE ( surf_usm_v(l)%r_a(1:surf_usm_v(l)%ns) ) 1303 ALLOCATE ( surf_usm_v(l)%r_a_green(1:surf_usm_v(l)%ns) ) 1304 ALLOCATE ( surf_usm_v(l)%r_a_window(1:surf_usm_v(l)%ns) ) 1277 1305 ALLOCATE ( surf_usm_v(l)%r_s(1:surf_usm_v(l)%ns) ) 1278 1306 ENDDO … … 1280 1308 ! 1281 1309 !-- allocate wall and roof material parameters. First for horizontal surfaces 1282 ALLOCATE ( surf_usm_h%thickness_wall(1:surf_usm_h%ns) )1310 ALLOCATE ( surf_usm_h%thickness_wall(1:surf_usm_h%ns) ) 1283 1311 ALLOCATE ( surf_usm_h%thickness_window(1:surf_usm_h%ns) ) 1284 1312 ALLOCATE ( surf_usm_h%thickness_green(1:surf_usm_h%ns) ) 1285 ALLOCATE ( surf_usm_h%lambda_h(nzb_wall:nzt_wall,1:surf_usm_h%ns) )1286 ALLOCATE ( surf_usm_h%rho_c_wall(nzb_wall:nzt_wall,1:surf_usm_h%ns) )1313 ALLOCATE ( surf_usm_h%lambda_h(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) 1314 ALLOCATE ( surf_usm_h%rho_c_wall(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) 1287 1315 ALLOCATE ( surf_usm_h%lambda_h_window(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) 1288 1316 ALLOCATE ( surf_usm_h%rho_c_window(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) … … 1290 1318 ALLOCATE ( surf_usm_h%rho_c_green(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) 1291 1319 1292 ALLOCATE ( surf_usm_h%rho_c_total_green(nzb_wall:nzt_wall,1:surf_usm_h%ns) 1293 ALLOCATE ( surf_usm_h%n_vg_green(1:surf_usm_h%ns) )1294 ALLOCATE ( surf_usm_h%alpha_vg_green(1:surf_usm_h%ns) )1295 ALLOCATE ( surf_usm_h%l_vg_green(1:surf_usm_h%ns) )1296 ALLOCATE ( surf_usm_h%gamma_w_green_sat(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) 1297 ALLOCATE ( surf_usm_h%lambda_w_green(nzb_wall:nzt_wall,1:surf_usm_h%ns) )1298 ALLOCATE ( surf_usm_h%gamma_w_green(nzb_wall:nzt_wall,1:surf_usm_h%ns) )1299 ALLOCATE ( surf_usm_h%tswc_h_m(nzb_wall:nzt_wall,1:surf_usm_h%ns) )1320 ALLOCATE ( surf_usm_h%rho_c_total_green(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) 1321 ALLOCATE ( surf_usm_h%n_vg_green(1:surf_usm_h%ns) ) 1322 ALLOCATE ( surf_usm_h%alpha_vg_green(1:surf_usm_h%ns) ) 1323 ALLOCATE ( surf_usm_h%l_vg_green(1:surf_usm_h%ns) ) 1324 ALLOCATE ( surf_usm_h%gamma_w_green_sat(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 1325 ALLOCATE ( surf_usm_h%lambda_w_green(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) 1326 ALLOCATE ( surf_usm_h%gamma_w_green(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) 1327 ALLOCATE ( surf_usm_h%tswc_h_m(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) 1300 1328 1301 1329 ! 1302 1330 !-- For vertical surfaces. 1303 1331 DO l = 0, 3 1304 ALLOCATE ( surf_usm_v(l)%thickness_wall(1:surf_usm_v(l)%ns) )1332 ALLOCATE ( surf_usm_v(l)%thickness_wall(1:surf_usm_v(l)%ns) ) 1305 1333 ALLOCATE ( surf_usm_v(l)%thickness_window(1:surf_usm_v(l)%ns) ) 1306 1334 ALLOCATE ( surf_usm_v(l)%thickness_green(1:surf_usm_v(l)%ns) ) 1307 ALLOCATE ( surf_usm_v(l)%lambda_h(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )1308 ALLOCATE ( surf_usm_v(l)%rho_c_wall(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )1335 ALLOCATE ( surf_usm_v(l)%lambda_h(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) ) 1336 ALLOCATE ( surf_usm_v(l)%rho_c_wall(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) ) 1309 1337 ALLOCATE ( surf_usm_v(l)%lambda_h_window(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) ) 1310 1338 ALLOCATE ( surf_usm_v(l)%rho_c_window(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) ) … … 1335 1363 ALLOCATE ( surf_usm_v(l)%r_canopy(1:surf_usm_v(l)%ns) ) 1336 1364 ALLOCATE ( surf_usm_v(l)%r_canopy_min(1:surf_usm_v(l)%ns) ) 1337 ALLOCATE ( surf_usm_v(l)%pt_10cm(1:surf_usm_v(l)%ns) )1365 ALLOCATE ( surf_usm_v(l)%pt_10cm(1:surf_usm_v(l)%ns) ) 1338 1366 ENDDO 1339 1367 1368 ! 1340 1369 !-- allocate wall and roof layers sizes. For horizontal surfaces. 1341 ALLOCATE ( zwn(nzb_wall:nzt_wall) )1342 ALLOCATE ( surf_usm_h%dz_wall(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )1343 ALLOCATE ( zwn_window(nzb_wall:nzt_wall) )1370 ALLOCATE ( zwn(nzb_wall:nzt_wall) ) 1371 ALLOCATE ( surf_usm_h%dz_wall(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 1372 ALLOCATE ( zwn_window(nzb_wall:nzt_wall) ) 1344 1373 ALLOCATE ( surf_usm_h%dz_window(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 1345 ALLOCATE ( zwn_green(nzb_wall:nzt_wall) )1374 ALLOCATE ( zwn_green(nzb_wall:nzt_wall) ) 1346 1375 ALLOCATE ( surf_usm_h%dz_green(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 1347 ALLOCATE ( surf_usm_h%ddz_wall(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )1348 ALLOCATE ( surf_usm_h%dz_wall_stag(nzb_wall:nzt_wall,1:surf_usm_h%ns) )1349 ALLOCATE ( surf_usm_h%ddz_wall_stag(nzb_wall:nzt_wall,1:surf_usm_h%ns) )1350 ALLOCATE ( surf_usm_h%zw(nzb_wall:nzt_wall,1:surf_usm_h%ns) )1376 ALLOCATE ( surf_usm_h%ddz_wall(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 1377 ALLOCATE ( surf_usm_h%dz_wall_stag(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) 1378 ALLOCATE ( surf_usm_h%ddz_wall_stag(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) 1379 ALLOCATE ( surf_usm_h%zw(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) 1351 1380 ALLOCATE ( surf_usm_h%ddz_window(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 1352 1381 ALLOCATE ( surf_usm_h%dz_window_stag(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) … … 1357 1386 ALLOCATE ( surf_usm_h%ddz_green_stag(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) 1358 1387 ALLOCATE ( surf_usm_h%zw_green(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) 1388 1359 1389 ! 1360 1390 !-- For vertical surfaces. 1361 1391 DO l = 0, 3 1362 ALLOCATE ( surf_usm_v(l)%dz_wall(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )1392 ALLOCATE ( surf_usm_v(l)%dz_wall(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) ) 1363 1393 ALLOCATE ( surf_usm_v(l)%dz_window(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) ) 1364 1394 ALLOCATE ( surf_usm_v(l)%dz_green(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) ) 1365 ALLOCATE ( surf_usm_v(l)%ddz_wall(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )1366 ALLOCATE ( surf_usm_v(l)%dz_wall_stag(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )1367 ALLOCATE ( surf_usm_v(l)%ddz_wall_stag(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )1368 ALLOCATE ( surf_usm_v(l)%zw(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )1395 ALLOCATE ( surf_usm_v(l)%ddz_wall(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) ) 1396 ALLOCATE ( surf_usm_v(l)%dz_wall_stag(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) ) 1397 ALLOCATE ( surf_usm_v(l)%ddz_wall_stag(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) ) 1398 ALLOCATE ( surf_usm_v(l)%zw(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) ) 1369 1399 ALLOCATE ( surf_usm_v(l)%ddz_window(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) ) 1370 1400 ALLOCATE ( surf_usm_v(l)%dz_window_stag(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) ) … … 1377 1407 ENDDO 1378 1408 1409 ! 1379 1410 !-- allocate wall and roof temperature arrays, for horizontal walls 1380 1411 ! 1381 1412 !-- Allocate if required. Note, in case of restarts, some of these arrays 1382 1413 !-- might be already allocated. 1383 IF ( .NOT. ALLOCATED( t_surf_wall_h_1 ) ) 1414 IF ( .NOT. ALLOCATED( t_surf_wall_h_1 ) ) & 1384 1415 ALLOCATE ( t_surf_wall_h_1(1:surf_usm_h%ns) ) 1385 IF ( .NOT. ALLOCATED( t_surf_wall_h_2 ) ) 1416 IF ( .NOT. ALLOCATED( t_surf_wall_h_2 ) ) & 1386 1417 ALLOCATE ( t_surf_wall_h_2(1:surf_usm_h%ns) ) 1387 1418 IF ( .NOT. ALLOCATED( t_wall_h_1 ) ) & … … 1405 1436 IF ( .NOT. ALLOCATED( t_green_h_2 ) ) & 1406 1437 ALLOCATE ( t_green_h_2(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 1407 IF ( .NOT. ALLOCATED( swc_h_1 ) ) &1438 IF ( .NOT. ALLOCATED( swc_h_1 ) ) & 1408 1439 ALLOCATE ( swc_h_1(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 1409 1440 IF ( .NOT. ALLOCATED( swc_sat_h_1 ) ) & … … 1411 1442 IF ( .NOT. ALLOCATED( swc_res_h_1 ) ) & 1412 1443 ALLOCATE ( swc_res_h_1(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 1413 IF ( .NOT. ALLOCATED( swc_h_2 ) ) &1444 IF ( .NOT. ALLOCATED( swc_h_2 ) ) & 1414 1445 ALLOCATE ( swc_h_2(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 1415 IF ( .NOT. ALLOCATED( rootfr_h_1 ) ) &1446 IF ( .NOT. ALLOCATED( rootfr_h_1 ) ) & 1416 1447 ALLOCATE ( rootfr_h_1(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 1417 IF ( .NOT. ALLOCATED( wilt_h_1 ) ) &1448 IF ( .NOT. ALLOCATED( wilt_h_1 ) ) & 1418 1449 ALLOCATE ( wilt_h_1(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 1419 IF ( .NOT. ALLOCATED( fc_h_1 ) ) &1450 IF ( .NOT. ALLOCATED( fc_h_1 ) ) & 1420 1451 ALLOCATE ( fc_h_1(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 1421 1452 1422 IF ( .NOT. ALLOCATED( m_liq_usm_h_1%var_usm_1d ) ) 1453 IF ( .NOT. ALLOCATED( m_liq_usm_h_1%var_usm_1d ) ) & 1423 1454 ALLOCATE ( m_liq_usm_h_1%var_usm_1d(1:surf_usm_h%ns) ) 1424 IF ( .NOT. ALLOCATED( m_liq_usm_h_2%var_usm_1d ) ) 1455 IF ( .NOT. ALLOCATED( m_liq_usm_h_2%var_usm_1d ) ) & 1425 1456 ALLOCATE ( m_liq_usm_h_2%var_usm_1d(1:surf_usm_h%ns) ) 1426 1457 1427 1458 ! 1428 1459 !-- initial assignment of the pointers 1429 t_wall_h => t_wall_h_1; t_wall_h_p=> t_wall_h_21430 t_window_h => t_window_h_1; t_window_h_p=> t_window_h_21431 t_green_h => t_green_h_1; t_green_h_p=> t_green_h_21432 t_surf_wall_h => t_surf_wall_h_1; t_surf_wall_h_p=> t_surf_wall_h_21460 t_wall_h => t_wall_h_1; t_wall_h_p => t_wall_h_2 1461 t_window_h => t_window_h_1; t_window_h_p => t_window_h_2 1462 t_green_h => t_green_h_1; t_green_h_p => t_green_h_2 1463 t_surf_wall_h => t_surf_wall_h_1; t_surf_wall_h_p => t_surf_wall_h_2 1433 1464 t_surf_window_h => t_surf_window_h_1; t_surf_window_h_p => t_surf_window_h_2 1434 t_surf_green_h => t_surf_green_h_1; t_surf_green_h_p=> t_surf_green_h_21465 t_surf_green_h => t_surf_green_h_1; t_surf_green_h_p => t_surf_green_h_2 1435 1466 m_liq_usm_h => m_liq_usm_h_1; m_liq_usm_h_p => m_liq_usm_h_2 1436 swc_h => swc_h_1; swc_h_p => swc_h_2 1437 swc_sat_h => swc_sat_h_1 1438 swc_res_h => swc_res_h_1 1439 rootfr_h => rootfr_h_1 1440 wilt_h => wilt_h_1 1441 fc_h => fc_h_1 1442 1467 swc_h => swc_h_1; swc_h_p => swc_h_2 1468 swc_sat_h => swc_sat_h_1 1469 swc_res_h => swc_res_h_1 1470 rootfr_h => rootfr_h_1 1471 wilt_h => wilt_h_1 1472 fc_h => fc_h_1 1473 1474 ! 1443 1475 !-- allocate wall and roof temperature arrays, for vertical walls if required 1444 1476 ! … … 1446 1478 !-- might be already allocated. 1447 1479 DO l = 0, 3 1448 IF ( .NOT. ALLOCATED( t_surf_wall_v_1(l)%t ) ) 1480 IF ( .NOT. ALLOCATED( t_surf_wall_v_1(l)%t ) ) & 1449 1481 ALLOCATE ( t_surf_wall_v_1(l)%t(1:surf_usm_v(l)%ns) ) 1450 IF ( .NOT. ALLOCATED( t_surf_wall_v_2(l)%t ) ) 1482 IF ( .NOT. ALLOCATED( t_surf_wall_v_2(l)%t ) ) & 1451 1483 ALLOCATE ( t_surf_wall_v_2(l)%t(1:surf_usm_v(l)%ns) ) 1452 1484 IF ( .NOT. ALLOCATED( t_wall_v_1(l)%t ) ) & … … 1470 1502 IF ( .NOT. ALLOCATED( t_green_v_2(l)%t ) ) & 1471 1503 ALLOCATE ( t_green_v_2(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) ) 1472 IF ( .NOT. ALLOCATED( m_liq_usm_v_1(l)%var_usm_1d ) ) &1504 IF ( .NOT. ALLOCATED( m_liq_usm_v_1(l)%var_usm_1d ) ) & 1473 1505 ALLOCATE ( m_liq_usm_v_1(l)%var_usm_1d(1:surf_usm_v(l)%ns) ) 1474 IF ( .NOT. ALLOCATED( m_liq_usm_v_2(l)%var_usm_1d ) ) &1506 IF ( .NOT. ALLOCATED( m_liq_usm_v_2(l)%var_usm_1d ) ) & 1475 1507 ALLOCATE ( m_liq_usm_v_2(l)%var_usm_1d(1:surf_usm_v(l)%ns) ) 1476 IF ( .NOT. ALLOCATED( swc_v_1(l)%t ) ) &1508 IF ( .NOT. ALLOCATED( swc_v_1(l)%t ) ) & 1477 1509 ALLOCATE ( swc_v_1(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) ) 1478 IF ( .NOT. ALLOCATED( swc_v_2(l)%t ) ) &1510 IF ( .NOT. ALLOCATED( swc_v_2(l)%t ) ) & 1479 1511 ALLOCATE ( swc_v_2(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) ) 1480 1512 ENDDO 1481 1513 ! 1482 1514 !-- initial assignment of the pointers 1483 t_wall_v => t_wall_v_1; t_wall_v_p=> t_wall_v_21484 t_surf_wall_v => t_surf_wall_v_1; t_surf_wall_v_p=> t_surf_wall_v_21485 t_window_v => t_window_v_1; t_window_v_p=> t_window_v_21486 t_green_v => t_green_v_1; t_green_v_p=> t_green_v_21515 t_wall_v => t_wall_v_1; t_wall_v_p => t_wall_v_2 1516 t_surf_wall_v => t_surf_wall_v_1; t_surf_wall_v_p => t_surf_wall_v_2 1517 t_window_v => t_window_v_1; t_window_v_p => t_window_v_2 1518 t_green_v => t_green_v_1; t_green_v_p => t_green_v_2 1487 1519 t_surf_window_v => t_surf_window_v_1; t_surf_window_v_p => t_surf_window_v_2 1488 t_surf_green_v => t_surf_green_v_1; t_surf_green_v_p=> t_surf_green_v_21520 t_surf_green_v => t_surf_green_v_1; t_surf_green_v_p => t_surf_green_v_2 1489 1521 m_liq_usm_v => m_liq_usm_v_1; m_liq_usm_v_p => m_liq_usm_v_2 1490 swc_v => swc_v_1; swc_v_p=> swc_v_21522 swc_v => swc_v_1; swc_v_p => swc_v_2 1491 1523 1492 1524 ! 1493 1525 !-- Allocate intermediate timestep arrays. For horizontal surfaces. 1494 ALLOCATE ( surf_usm_h%tt_surface_wall_m(1:surf_usm_h%ns) 1495 ALLOCATE ( surf_usm_h%tt_wall_m(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )1526 ALLOCATE ( surf_usm_h%tt_surface_wall_m(1:surf_usm_h%ns) ) 1527 ALLOCATE ( surf_usm_h%tt_wall_m(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 1496 1528 ALLOCATE ( surf_usm_h%tt_surface_window_m(1:surf_usm_h%ns) ) 1497 1529 ALLOCATE ( surf_usm_h%tt_window_m(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) … … 1511 1543 ! 1512 1544 !-- Set inital values for prognostic quantities 1513 IF ( ALLOCATED( surf_usm_h%tt_surface_wall_m ) ) surf_usm_h%tt_surface_wall_m= 0.0_wp1514 IF ( ALLOCATED( surf_usm_h%tt_wall_m ) ) surf_usm_h%tt_wall_m= 0.0_wp1545 IF ( ALLOCATED( surf_usm_h%tt_surface_wall_m ) ) surf_usm_h%tt_surface_wall_m = 0.0_wp 1546 IF ( ALLOCATED( surf_usm_h%tt_wall_m ) ) surf_usm_h%tt_wall_m = 0.0_wp 1515 1547 IF ( ALLOCATED( surf_usm_h%tt_surface_window_m ) ) surf_usm_h%tt_surface_window_m = 0.0_wp 1516 1548 IF ( ALLOCATED( surf_usm_h%tt_window_m ) ) surf_usm_h%tt_window_m = 0.0_wp … … 1520 1552 !-- Now, for vertical surfaces 1521 1553 DO l = 0, 3 1522 ALLOCATE ( surf_usm_v(l)%tt_surface_wall_m(1:surf_usm_v(l)%ns) 1523 ALLOCATE ( surf_usm_v(l)%tt_wall_m(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )1554 ALLOCATE ( surf_usm_v(l)%tt_surface_wall_m(1:surf_usm_v(l)%ns) ) 1555 ALLOCATE ( surf_usm_v(l)%tt_wall_m(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) ) 1524 1556 IF ( ALLOCATED( surf_usm_v(l)%tt_surface_wall_m ) ) surf_usm_v(l)%tt_surface_wall_m = 0.0_wp 1525 1557 IF ( ALLOCATED( surf_usm_v(l)%tt_wall_m ) ) surf_usm_v(l)%tt_wall_m = 0.0_wp … … 1533 1565 IF ( ALLOCATED( surf_usm_v(l)%tt_green_m ) ) surf_usm_v(l)%tt_green_m = 0.0_wp 1534 1566 ENDDO 1535 1567 ! 1536 1568 !-- allocate wall heat flux output array and set initial values. For horizontal surfaces 1537 ! 1569 ! ALLOCATE ( surf_usm_h%wshf(1:surf_usm_h%ns) ) !can be removed 1538 1570 ALLOCATE ( surf_usm_h%wshf_eb(1:surf_usm_h%ns) ) 1539 1571 ALLOCATE ( surf_usm_h%wghf_eb(1:surf_usm_h%ns) ) … … 1552 1584 !-- Now, for vertical surfaces 1553 1585 DO l = 0, 3 1554 ! 1586 ! ALLOCATE ( surf_usm_v(l)%wshf(1:surf_usm_v(l)%ns) ) ! can be removed 1555 1587 ALLOCATE ( surf_usm_v(l)%wshf_eb(1:surf_usm_v(l)%ns) ) 1556 1588 ALLOCATE ( surf_usm_v(l)%wghf_eb(1:surf_usm_v(l)%ns) ) … … 1569 1601 1570 1602 CALL location_message( 'finished', .TRUE. ) 1571 1603 1572 1604 END SUBROUTINE usm_init_arrays 1573 1605 … … 1586 1618 CHARACTER(LEN=*), INTENT(IN) :: variable 1587 1619 1588 INTEGER(iwp) :: i, j, k, l, m, ids, idsint, iwl, istat 1589 CHARACTER(LEN=varnamelength) :: var 1590 INTEGER(iwp), PARAMETER :: nd = 5 1620 INTEGER(iwp) :: i, j, k, l, m, ids, idsint, iwl, istat !< runnin indices 1621 CHARACTER(LEN=varnamelength) :: var !< trimmed variable 1622 INTEGER(iwp), PARAMETER :: nd = 5 !< number of directions 1591 1623 CHARACTER(LEN=6), DIMENSION(0:nd-1), PARAMETER :: dirname = (/ '_roof ', '_south', '_north', '_west ', '_east ' /) 1592 1624 INTEGER(iwp), DIMENSION(0:nd-1), PARAMETER :: dirint = (/ iup_u, isouth_u, inorth_u, iwest_u, ieast_u /) … … 1594 1626 IF ( variable(1:4) == 'usm_' ) THEN ! is such a check really rquired? 1595 1627 1628 ! 1596 1629 !-- find the real name of the variable 1597 1630 ids = -1 … … 1616 1649 ENDIF 1617 1650 IF ( var(1:11) == 'usm_t_wall_' .AND. len(TRIM(var)) >= 12 ) THEN 1651 ! 1618 1652 !-- wall layers 1619 1653 READ(var(12:12), '(I1)', iostat=istat ) iwl … … 1621 1655 var = var(1:10) 1622 1656 ELSE 1657 ! 1623 1658 !-- wrong wall layer index 1624 1659 RETURN … … 1626 1661 ENDIF 1627 1662 IF ( var(1:13) == 'usm_t_window_' .AND. len(TRIM(var)) >= 14 ) THEN 1663 ! 1628 1664 !-- wall layers 1629 1665 READ(var(14:14), '(I1)', iostat=istat ) iwl … … 1631 1667 var = var(1:12) 1632 1668 ELSE 1669 ! 1633 1670 !-- wrong window layer index 1634 1671 RETURN … … 1636 1673 ENDIF 1637 1674 IF ( var(1:12) == 'usm_t_green_' .AND. len(TRIM(var)) >= 13 ) THEN 1675 ! 1638 1676 !-- wall layers 1639 1677 READ(var(13:13), '(I1)', iostat=istat ) iwl … … 1641 1679 var = var(1:11) 1642 1680 ELSE 1681 ! 1643 1682 !-- wrong green layer index 1644 1683 RETURN … … 1646 1685 ENDIF 1647 1686 IF ( var(1:8) == 'usm_swc_' .AND. len(TRIM(var)) >= 9 ) THEN 1687 ! 1648 1688 !-- swc layers 1649 1689 READ(var(9:9), '(I1)', iostat=istat ) iwl … … 1651 1691 var = var(1:7) 1652 1692 ELSE 1693 ! 1653 1694 !-- wrong swc layer index 1654 1695 RETURN … … 1661 1702 1662 1703 CASE ( 'usm_wshf' ) 1704 ! 1663 1705 !-- array of sensible heat flux from surfaces 1664 1706 !-- land surfaces 1665 1707 IF ( l == -1 ) THEN 1666 1708 IF ( .NOT. ALLOCATED(surf_usm_h%wshf_eb_av) ) THEN 1667 ALLOCATE ( surf_usm_h%wshf_eb_av(1:surf_usm_h%ns) )1709 ALLOCATE ( surf_usm_h%wshf_eb_av(1:surf_usm_h%ns) ) 1668 1710 surf_usm_h%wshf_eb_av = 0.0_wp 1669 1711 ENDIF 1670 1712 ELSE 1671 1713 IF ( .NOT. ALLOCATED(surf_usm_v(l)%wshf_eb_av) ) THEN 1672 ALLOCATE ( surf_usm_v(l)%wshf_eb_av(1:surf_usm_v(l)%ns) )1714 ALLOCATE ( surf_usm_v(l)%wshf_eb_av(1:surf_usm_v(l)%ns) ) 1673 1715 surf_usm_v(l)%wshf_eb_av = 0.0_wp 1674 1716 ENDIF … … 1676 1718 1677 1719 CASE ( 'usm_qsws' ) 1720 ! 1678 1721 !-- array of latent heat flux from surfaces 1679 1722 !-- land surfaces 1680 1723 IF ( l == -1 .AND. .NOT. ALLOCATED(surf_usm_h%qsws_eb_av) ) THEN 1681 ALLOCATE ( surf_usm_h%qsws_eb_av(1:surf_usm_h%ns) )1724 ALLOCATE ( surf_usm_h%qsws_eb_av(1:surf_usm_h%ns) ) 1682 1725 surf_usm_h%qsws_eb_av = 0.0_wp 1683 1726 ELSE 1684 1727 IF ( .NOT. ALLOCATED(surf_usm_v(l)%qsws_eb_av) ) THEN 1685 ALLOCATE ( surf_usm_v(l)%qsws_eb_av(1:surf_usm_v(l)%ns) )1728 ALLOCATE ( surf_usm_v(l)%qsws_eb_av(1:surf_usm_v(l)%ns) ) 1686 1729 surf_usm_v(l)%qsws_eb_av = 0.0_wp 1687 1730 ENDIF … … 1689 1732 1690 1733 CASE ( 'usm_qsws_veg' ) 1734 ! 1691 1735 !-- array of latent heat flux from vegetation surfaces 1692 1736 !-- land surfaces 1693 1737 IF ( l == -1 .AND. .NOT. ALLOCATED(surf_usm_h%qsws_veg_av) ) THEN 1694 ALLOCATE ( surf_usm_h%qsws_veg_av(1:surf_usm_h%ns) )1738 ALLOCATE ( surf_usm_h%qsws_veg_av(1:surf_usm_h%ns) ) 1695 1739 surf_usm_h%qsws_veg_av = 0.0_wp 1696 1740 ELSE 1697 1741 IF ( .NOT. ALLOCATED(surf_usm_v(l)%qsws_veg_av) ) THEN 1698 ALLOCATE ( surf_usm_v(l)%qsws_veg_av(1:surf_usm_v(l)%ns) )1742 ALLOCATE ( surf_usm_v(l)%qsws_veg_av(1:surf_usm_v(l)%ns) ) 1699 1743 surf_usm_v(l)%qsws_veg_av = 0.0_wp 1700 1744 ENDIF … … 1702 1746 1703 1747 CASE ( 'usm_qsws_liq' ) 1748 ! 1704 1749 !-- array of latent heat flux from surfaces with liquid 1705 1750 !-- land surfaces 1706 1751 IF ( l == -1 .AND. .NOT. ALLOCATED(surf_usm_h%qsws_liq_av) ) THEN 1707 ALLOCATE ( surf_usm_h%qsws_liq_av(1:surf_usm_h%ns) )1752 ALLOCATE ( surf_usm_h%qsws_liq_av(1:surf_usm_h%ns) ) 1708 1753 surf_usm_h%qsws_liq_av = 0.0_wp 1709 1754 ELSE 1710 1755 IF ( .NOT. ALLOCATED(surf_usm_v(l)%qsws_liq_av) ) THEN 1711 ALLOCATE ( surf_usm_v(l)%qsws_liq_av(1:surf_usm_v(l)%ns) )1756 ALLOCATE ( surf_usm_v(l)%qsws_liq_av(1:surf_usm_v(l)%ns) ) 1712 1757 surf_usm_v(l)%qsws_liq_av = 0.0_wp 1713 1758 ENDIF … … 1719 1764 !-- accordingly in average_3d_data, sum_up_3d_data, etc.. 1720 1765 CASE ( 'usm_wghf' ) 1766 ! 1721 1767 !-- array of heat flux from ground (wall, roof, land) 1722 1768 IF ( l == -1 ) THEN 1723 1769 IF ( .NOT. ALLOCATED(surf_usm_h%wghf_eb_av) ) THEN 1724 ALLOCATE ( surf_usm_h%wghf_eb_av(1:surf_usm_h%ns) )1770 ALLOCATE ( surf_usm_h%wghf_eb_av(1:surf_usm_h%ns) ) 1725 1771 surf_usm_h%wghf_eb_av = 0.0_wp 1726 1772 ENDIF 1727 1773 ELSE 1728 1774 IF ( .NOT. ALLOCATED(surf_usm_v(l)%wghf_eb_av) ) THEN 1729 ALLOCATE ( surf_usm_v(l)%wghf_eb_av(1:surf_usm_v(l)%ns) )1775 ALLOCATE ( surf_usm_v(l)%wghf_eb_av(1:surf_usm_v(l)%ns) ) 1730 1776 surf_usm_v(l)%wghf_eb_av = 0.0_wp 1731 1777 ENDIF … … 1733 1779 1734 1780 CASE ( 'usm_wghf_window' ) 1781 ! 1735 1782 !-- array of heat flux from window ground (wall, roof, land) 1736 1783 IF ( l == -1 ) THEN 1737 1784 IF ( .NOT. ALLOCATED(surf_usm_h%wghf_eb_window_av) ) THEN 1738 ALLOCATE ( surf_usm_h%wghf_eb_window_av(1:surf_usm_h%ns) )1785 ALLOCATE ( surf_usm_h%wghf_eb_window_av(1:surf_usm_h%ns) ) 1739 1786 surf_usm_h%wghf_eb_window_av = 0.0_wp 1740 1787 ENDIF 1741 1788 ELSE 1742 1789 IF ( .NOT. ALLOCATED(surf_usm_v(l)%wghf_eb_window_av) ) THEN 1743 ALLOCATE ( surf_usm_v(l)%wghf_eb_window_av(1:surf_usm_v(l)%ns) )1790 ALLOCATE ( surf_usm_v(l)%wghf_eb_window_av(1:surf_usm_v(l)%ns) ) 1744 1791 surf_usm_v(l)%wghf_eb_window_av = 0.0_wp 1745 1792 ENDIF … … 1747 1794 1748 1795 CASE ( 'usm_wghf_green' ) 1796 ! 1749 1797 !-- array of heat flux from green ground (wall, roof, land) 1750 1798 IF ( l == -1 ) THEN 1751 1799 IF ( .NOT. ALLOCATED(surf_usm_h%wghf_eb_green_av) ) THEN 1752 ALLOCATE ( surf_usm_h%wghf_eb_green_av(1:surf_usm_h%ns) )1800 ALLOCATE ( surf_usm_h%wghf_eb_green_av(1:surf_usm_h%ns) ) 1753 1801 surf_usm_h%wghf_eb_green_av = 0.0_wp 1754 1802 ENDIF 1755 1803 ELSE 1756 1804 IF ( .NOT. ALLOCATED(surf_usm_v(l)%wghf_eb_green_av) ) THEN 1757 ALLOCATE ( surf_usm_v(l)%wghf_eb_green_av(1:surf_usm_v(l)%ns) )1805 ALLOCATE ( surf_usm_v(l)%wghf_eb_green_av(1:surf_usm_v(l)%ns) ) 1758 1806 surf_usm_v(l)%wghf_eb_green_av = 0.0_wp 1759 1807 ENDIF … … 1761 1809 1762 1810 CASE ( 'usm_iwghf' ) 1811 ! 1763 1812 !-- array of heat flux from indoor ground (wall, roof, land) 1764 1813 IF ( l == -1 ) THEN 1765 1814 IF ( .NOT. ALLOCATED(surf_usm_h%iwghf_eb_av) ) THEN 1766 ALLOCATE ( surf_usm_h%iwghf_eb_av(1:surf_usm_h%ns) )1815 ALLOCATE ( surf_usm_h%iwghf_eb_av(1:surf_usm_h%ns) ) 1767 1816 surf_usm_h%iwghf_eb_av = 0.0_wp 1768 1817 ENDIF 1769 1818 ELSE 1770 1819 IF ( .NOT. ALLOCATED(surf_usm_v(l)%iwghf_eb_av) ) THEN 1771 ALLOCATE ( surf_usm_v(l)%iwghf_eb_av(1:surf_usm_v(l)%ns) )1820 ALLOCATE ( surf_usm_v(l)%iwghf_eb_av(1:surf_usm_v(l)%ns) ) 1772 1821 surf_usm_v(l)%iwghf_eb_av = 0.0_wp 1773 1822 ENDIF … … 1775 1824 1776 1825 CASE ( 'usm_iwghf_window' ) 1826 ! 1777 1827 !-- array of heat flux from indoor window ground (wall, roof, land) 1778 1828 IF ( l == -1 ) THEN 1779 1829 IF ( .NOT. ALLOCATED(surf_usm_h%iwghf_eb_window_av) ) THEN 1780 ALLOCATE ( surf_usm_h%iwghf_eb_window_av(1:surf_usm_h%ns) )1830 ALLOCATE ( surf_usm_h%iwghf_eb_window_av(1:surf_usm_h%ns) ) 1781 1831 surf_usm_h%iwghf_eb_window_av = 0.0_wp 1782 1832 ENDIF 1783 1833 ELSE 1784 1834 IF ( .NOT. ALLOCATED(surf_usm_v(l)%iwghf_eb_window_av) ) THEN 1785 ALLOCATE ( surf_usm_v(l)%iwghf_eb_window_av(1:surf_usm_v(l)%ns) )1835 ALLOCATE ( surf_usm_v(l)%iwghf_eb_window_av(1:surf_usm_v(l)%ns) ) 1786 1836 surf_usm_v(l)%iwghf_eb_window_av = 0.0_wp 1787 1837 ENDIF … … 1789 1839 1790 1840 CASE ( 'usm_t_surf_wall' ) 1841 ! 1791 1842 !-- surface temperature for surfaces 1792 1843 IF ( l == -1 ) THEN 1793 1844 IF ( .NOT. ALLOCATED(surf_usm_h%t_surf_wall_av) ) THEN 1794 ALLOCATE ( surf_usm_h%t_surf_wall_av(1:surf_usm_h%ns) )1845 ALLOCATE ( surf_usm_h%t_surf_wall_av(1:surf_usm_h%ns) ) 1795 1846 surf_usm_h%t_surf_wall_av = 0.0_wp 1796 1847 ENDIF 1797 1848 ELSE 1798 1849 IF ( .NOT. ALLOCATED(surf_usm_v(l)%t_surf_wall_av) ) THEN 1799 ALLOCATE ( surf_usm_v(l)%t_surf_wall_av(1:surf_usm_v(l)%ns) )1850 ALLOCATE ( surf_usm_v(l)%t_surf_wall_av(1:surf_usm_v(l)%ns) ) 1800 1851 surf_usm_v(l)%t_surf_wall_av = 0.0_wp 1801 1852 ENDIF … … 1803 1854 1804 1855 CASE ( 'usm_t_surf_window' ) 1856 ! 1805 1857 !-- surface temperature for window surfaces 1806 1858 IF ( l == -1 ) THEN 1807 1859 IF ( .NOT. ALLOCATED(surf_usm_h%t_surf_window_av) ) THEN 1808 ALLOCATE ( surf_usm_h%t_surf_window_av(1:surf_usm_h%ns) )1860 ALLOCATE ( surf_usm_h%t_surf_window_av(1:surf_usm_h%ns) ) 1809 1861 surf_usm_h%t_surf_window_av = 0.0_wp 1810 1862 ENDIF 1811 1863 ELSE 1812 1864 IF ( .NOT. ALLOCATED(surf_usm_v(l)%t_surf_window_av) ) THEN 1813 ALLOCATE ( surf_usm_v(l)%t_surf_window_av(1:surf_usm_v(l)%ns) )1865 ALLOCATE ( surf_usm_v(l)%t_surf_window_av(1:surf_usm_v(l)%ns) ) 1814 1866 surf_usm_v(l)%t_surf_window_av = 0.0_wp 1815 1867 ENDIF … … 1817 1869 1818 1870 CASE ( 'usm_t_surf_green' ) 1871 ! 1819 1872 !-- surface temperature for green surfaces 1820 1873 IF ( l == -1 ) THEN 1821 1874 IF ( .NOT. ALLOCATED(surf_usm_h%t_surf_green_av) ) THEN 1822 ALLOCATE ( surf_usm_h%t_surf_green_av(1:surf_usm_h%ns) )1875 ALLOCATE ( surf_usm_h%t_surf_green_av(1:surf_usm_h%ns) ) 1823 1876 surf_usm_h%t_surf_green_av = 0.0_wp 1824 1877 ENDIF 1825 1878 ELSE 1826 1879 IF ( .NOT. ALLOCATED(surf_usm_v(l)%t_surf_green_av) ) THEN 1827 ALLOCATE ( surf_usm_v(l)%t_surf_green_av(1:surf_usm_v(l)%ns) )1880 ALLOCATE ( surf_usm_v(l)%t_surf_green_av(1:surf_usm_v(l)%ns) ) 1828 1881 surf_usm_v(l)%t_surf_green_av = 0.0_wp 1829 1882 ENDIF … … 1831 1884 1832 1885 CASE ( 'usm_theta_10cm' ) 1886 ! 1833 1887 !-- near surface (10cm) temperature for whole surfaces 1834 1888 IF ( l == -1 ) THEN 1835 1889 IF ( .NOT. ALLOCATED(surf_usm_h%pt_10cm_av) ) THEN 1836 ALLOCATE ( surf_usm_h%pt_10cm_av(1:surf_usm_h%ns) )1890 ALLOCATE ( surf_usm_h%pt_10cm_av(1:surf_usm_h%ns) ) 1837 1891 surf_usm_h%pt_10cm_av = 0.0_wp 1838 1892 ENDIF 1839 1893 ELSE 1840 1894 IF ( .NOT. ALLOCATED(surf_usm_v(l)%pt_10cm_av) ) THEN 1841 ALLOCATE ( surf_usm_v(l)%pt_10cm_av(1:surf_usm_v(l)%ns) )1895 ALLOCATE ( surf_usm_v(l)%pt_10cm_av(1:surf_usm_v(l)%ns) ) 1842 1896 surf_usm_v(l)%pt_10cm_av = 0.0_wp 1843 1897 ENDIF … … 1845 1899 1846 1900 CASE ( 'usm_t_wall' ) 1901 ! 1847 1902 !-- wall temperature for iwl layer of walls and land 1848 1903 IF ( l == -1 ) THEN 1849 1904 IF ( .NOT. ALLOCATED(surf_usm_h%t_wall_av) ) THEN 1850 ALLOCATE ( surf_usm_h%t_wall_av(nzb_wall:nzt_wall,1:surf_usm_h%ns) )1905 ALLOCATE ( surf_usm_h%t_wall_av(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) 1851 1906 surf_usm_h%t_wall_av = 0.0_wp 1852 1907 ENDIF 1853 1908 ELSE 1854 1909 IF ( .NOT. ALLOCATED(surf_usm_v(l)%t_wall_av) ) THEN 1855 ALLOCATE ( surf_usm_v(l)%t_wall_av(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )1910 ALLOCATE ( surf_usm_v(l)%t_wall_av(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) ) 1856 1911 surf_usm_v(l)%t_wall_av = 0.0_wp 1857 1912 ENDIF … … 1859 1914 1860 1915 CASE ( 'usm_t_window' ) 1916 ! 1861 1917 !-- window temperature for iwl layer of walls and land 1862 1918 IF ( l == -1 ) THEN 1863 1919 IF ( .NOT. ALLOCATED(surf_usm_h%t_window_av) ) THEN 1864 ALLOCATE ( surf_usm_h%t_window_av(nzb_wall:nzt_wall,1:surf_usm_h%ns) )1920 ALLOCATE ( surf_usm_h%t_window_av(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) 1865 1921 surf_usm_h%t_window_av = 0.0_wp 1866 1922 ENDIF 1867 1923 ELSE 1868 1924 IF ( .NOT. ALLOCATED(surf_usm_v(l)%t_window_av) ) THEN 1869 ALLOCATE ( surf_usm_v(l)%t_window_av(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )1925 ALLOCATE ( surf_usm_v(l)%t_window_av(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) ) 1870 1926 surf_usm_v(l)%t_window_av = 0.0_wp 1871 1927 ENDIF … … 1873 1929 1874 1930 CASE ( 'usm_t_green' ) 1931 ! 1875 1932 !-- green temperature for iwl layer of walls and land 1876 1933 IF ( l == -1 ) THEN 1877 1934 IF ( .NOT. ALLOCATED(surf_usm_h%t_green_av) ) THEN 1878 ALLOCATE ( surf_usm_h%t_green_av(nzb_wall:nzt_wall,1:surf_usm_h%ns) )1935 ALLOCATE ( surf_usm_h%t_green_av(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) 1879 1936 surf_usm_h%t_green_av = 0.0_wp 1880 1937 ENDIF 1881 1938 ELSE 1882 1939 IF ( .NOT. ALLOCATED(surf_usm_v(l)%t_green_av) ) THEN 1883 ALLOCATE ( surf_usm_v(l)%t_green_av(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )1940 ALLOCATE ( surf_usm_v(l)%t_green_av(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) ) 1884 1941 surf_usm_v(l)%t_green_av = 0.0_wp 1885 1942 ENDIF 1886 1943 ENDIF 1887 1944 CASE ( 'usm_swc' ) 1945 ! 1888 1946 !-- soil water content for iwl layer of walls and land 1889 1947 IF ( l == -1 .AND. .NOT. ALLOCATED(surf_usm_h%swc_av) ) THEN 1890 ALLOCATE ( surf_usm_h%swc_av(nzb_wall:nzt_wall,1:surf_usm_h%ns) )1948 ALLOCATE ( surf_usm_h%swc_av(nzb_wall:nzt_wall,1:surf_usm_h%ns) ) 1891 1949 surf_usm_h%swc_av = 0.0_wp 1892 1950 ELSE 1893 1951 IF ( .NOT. ALLOCATED(surf_usm_v(l)%swc_av) ) THEN 1894 ALLOCATE ( surf_usm_v(l)%swc_av(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )1952 ALLOCATE ( surf_usm_v(l)%swc_av(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) ) 1895 1953 surf_usm_v(l)%swc_av = 0.0_wp 1896 1954 ENDIF … … 1907 1965 1908 1966 CASE ( 'usm_wshf' ) 1967 ! 1909 1968 !-- array of sensible heat flux from surfaces (land, roof, wall) 1910 1969 IF ( l == -1 ) THEN … … 1923 1982 1924 1983 CASE ( 'usm_qsws' ) 1984 ! 1925 1985 !-- array of latent heat flux from surfaces (land, roof, wall) 1926 1986 IF ( l == -1 ) THEN … … 1939 1999 1940 2000 CASE ( 'usm_qsws_veg' ) 2001 ! 1941 2002 !-- array of latent heat flux from vegetation surfaces (land, roof, wall) 1942 2003 IF ( l == -1 ) THEN … … 1955 2016 1956 2017 CASE ( 'usm_qsws_liq' ) 2018 ! 1957 2019 !-- array of latent heat flux from surfaces with liquid (land, roof, wall) 1958 2020 IF ( l == -1 ) THEN … … 1971 2033 1972 2034 CASE ( 'usm_wghf' ) 2035 ! 1973 2036 !-- array of heat flux from ground (wall, roof, land) 1974 2037 IF ( l == -1 ) THEN … … 1987 2050 1988 2051 CASE ( 'usm_wghf_window' ) 2052 ! 1989 2053 !-- array of heat flux from window ground (wall, roof, land) 1990 2054 IF ( l == -1 ) THEN … … 2003 2067 2004 2068 CASE ( 'usm_wghf_green' ) 2069 ! 2005 2070 !-- array of heat flux from green ground (wall, roof, land) 2006 2071 IF ( l == -1 ) THEN … … 2019 2084 2020 2085 CASE ( 'usm_iwghf' ) 2086 ! 2021 2087 !-- array of heat flux from indoor ground (wall, roof, land) 2022 2088 IF ( l == -1 ) THEN … … 2035 2101 2036 2102 CASE ( 'usm_iwghf_window' ) 2103 ! 2037 2104 !-- array of heat flux from indoor window ground (wall, roof, land) 2038 2105 IF ( l == -1 ) THEN … … 2051 2118 2052 2119 CASE ( 'usm_t_surf_wall' ) 2120 ! 2053 2121 !-- surface temperature for surfaces 2054 2122 IF ( l == -1 ) THEN … … 2067 2135 2068 2136 CASE ( 'usm_t_surf_window' ) 2137 ! 2069 2138 !-- surface temperature for window surfaces 2070 2139 IF ( l == -1 ) THEN … … 2083 2152 2084 2153 CASE ( 'usm_t_surf_green' ) 2154 ! 2085 2155 !-- surface temperature for green surfaces 2086 2156 IF ( l == -1 ) THEN … … 2099 2169 2100 2170 CASE ( 'usm_theta_10cm' ) 2171 ! 2101 2172 !-- near surface temperature for whole surfaces 2102 2173 IF ( l == -1 ) THEN … … 2115 2186 2116 2187 CASE ( 'usm_t_wall' ) 2188 ! 2117 2189 !-- wall temperature for iwl layer of walls and land 2118 2190 IF ( l == -1 ) THEN … … 2131 2203 2132 2204 CASE ( 'usm_t_window' ) 2205 ! 2133 2206 !-- window temperature for iwl layer of walls and land 2134 2207 IF ( l == -1 ) THEN … … 2147 2220 2148 2221 CASE ( 'usm_t_green' ) 2222 ! 2149 2223 !-- green temperature for iwl layer of walls and land 2150 2224 IF ( l == -1 ) THEN … … 2163 2237 2164 2238 CASE ( 'usm_swc' ) 2239 ! 2165 2240 !-- soil water content for iwl layer of walls and land 2166 2241 IF ( l == -1 ) THEN … … 2188 2263 2189 2264 CASE ( 'usm_wshf' ) 2265 ! 2190 2266 !-- array of sensible heat flux from surfaces (land, roof, wall) 2191 2267 IF ( l == -1 ) THEN … … 2204 2280 2205 2281 CASE ( 'usm_qsws' ) 2282 ! 2206 2283 !-- array of latent heat flux from surfaces (land, roof, wall) 2207 2284 IF ( l == -1 ) THEN … … 2220 2297 2221 2298 CASE ( 'usm_qsws_veg' ) 2299 ! 2222 2300 !-- array of latent heat flux from vegetation surfaces (land, roof, wall) 2223 2301 IF ( l == -1 ) THEN … … 2236 2314 2237 2315 CASE ( 'usm_qsws_liq' ) 2316 ! 2238 2317 !-- array of latent heat flux from surfaces with liquid (land, roof, wall) 2239 2318 IF ( l == -1 ) THEN … … 2252 2331 2253 2332 CASE ( 'usm_wghf' ) 2333 ! 2254 2334 !-- array of heat flux from ground (wall, roof, land) 2255 2335 IF ( l == -1 ) THEN … … 2268 2348 2269 2349 CASE ( 'usm_wghf_window' ) 2350 ! 2270 2351 !-- array of heat flux from window ground (wall, roof, land) 2271 2352 IF ( l == -1 ) THEN … … 2284 2365 2285 2366 CASE ( 'usm_wghf_green' ) 2367 ! 2286 2368 !-- array of heat flux from green ground (wall, roof, land) 2287 2369 IF ( l == -1 ) THEN … … 2300 2382 2301 2383 CASE ( 'usm_iwghf' ) 2384 ! 2302 2385 !-- array of heat flux from indoor ground (wall, roof, land) 2303 2386 IF ( l == -1 ) THEN … … 2316 2399 2317 2400 CASE ( 'usm_iwghf_window' ) 2401 ! 2318 2402 !-- array of heat flux from indoor window ground (wall, roof, land) 2319 2403 IF ( l == -1 ) THEN … … 2332 2416 2333 2417 CASE ( 'usm_t_surf_wall' ) 2418 ! 2334 2419 !-- surface temperature for surfaces 2335 2420 IF ( l == -1 ) THEN … … 2348 2433 2349 2434 CASE ( 'usm_t_surf_window' ) 2435 ! 2350 2436 !-- surface temperature for window surfaces 2351 2437 IF ( l == -1 ) THEN … … 2364 2450 2365 2451 CASE ( 'usm_t_surf_green' ) 2452 ! 2366 2453 !-- surface temperature for green surfaces 2367 2454 IF ( l == -1 ) THEN … … 2380 2467 2381 2468 CASE ( 'usm_theta_10cm' ) 2469 ! 2382 2470 !-- near surface temperature for whole surfaces 2383 2471 IF ( l == -1 ) THEN … … 2397 2485 2398 2486 CASE ( 'usm_t_wall' ) 2487 ! 2399 2488 !-- wall temperature for iwl layer of walls and land 2400 2489 IF ( l == -1 ) THEN … … 2413 2502 2414 2503 CASE ( 'usm_t_window' ) 2504 ! 2415 2505 !-- window temperature for iwl layer of walls and land 2416 2506 IF ( l == -1 ) THEN … … 2429 2519 2430 2520 CASE ( 'usm_t_green' ) 2521 ! 2431 2522 !-- green temperature for iwl layer of walls and land 2432 2523 IF ( l == -1 ) THEN … … 2445 2536 2446 2537 CASE ( 'usm_swc' ) 2538 ! 2447 2539 !-- soil water content for iwl layer of walls and land 2448 2540 IF ( l == -1 ) THEN … … 2526 2618 CHARACTER(LEN=*),INTENT(OUT) :: unit !< 2527 2619 2528 INTEGER(iwp) :: i,j,l !< index2620 INTEGER(iwp) :: i,j,l !< index 2529 2621 CHARACTER(LEN=2) :: ls 2530 CHARACTER(LEN=varnamelength) :: var !< TRIM(variable)2531 INTEGER(iwp), PARAMETER :: nl1 = 16 !< number of directional usm variables2532 CHARACTER(LEN=varnamelength), DIMENSION(nl1) :: varlist1 = & !< list of directional usm variables2622 CHARACTER(LEN=varnamelength) :: var !< TRIM(variable) 2623 INTEGER(iwp), PARAMETER :: nl1 = 16 !< number of directional usm variables 2624 CHARACTER(LEN=varnamelength), DIMENSION(nl1) :: varlist1 = & !< list of directional usm variables 2533 2625 (/'usm_wshf ', & 2534 2626 'usm_wghf ', & … … 2548 2640 'usm_theta_10cm '/) 2549 2641 2550 INTEGER(iwp), PARAMETER :: nl2 = 3 !< number of directional layer usm variables2551 CHARACTER(LEN=varnamelength), DIMENSION(nl2) :: varlist2 = & !< list of directional layer usm variables2642 INTEGER(iwp), PARAMETER :: nl2 = 3 !< number of directional layer usm variables 2643 CHARACTER(LEN=varnamelength), DIMENSION(nl2) :: varlist2 = & !< list of directional layer usm variables 2552 2644 (/'usm_t_wall ', & 2553 2645 'usm_t_window ', & … … 2564 2656 var = TRIM(variable) 2565 2657 2658 ! 2566 2659 !-- check if variable exists 2567 ! 2660 !-- directional variables 2568 2661 DO i = 1, nl1 2569 2662 DO j = 1, nd … … 2576 2669 ENDDO 2577 2670 IF ( lfound ) GOTO 10 2578 ! directional layer variables 2671 ! 2672 !-- directional layer variables 2579 2673 DO i = 1, nl2 2580 2674 DO j = 1, nd … … 2624 2718 !------------------------------------------------------------------------------! 2625 2719 SUBROUTINE usm_check_parameters 2626 2720 2627 2721 USE control_parameters, & 2628 2722 ONLY: bc_pt_b, bc_q_b, constant_flux_layer, large_scale_forcing, & 2629 2723 lsf_surf, topography 2630 2724 2631 2725 USE netcdf_data_input_mod, & 2632 2726 ONLY: building_type_f 2633 2727 2634 IMPLICIT NONE 2635 2728 IMPLICIT NONE 2729 2636 2730 INTEGER(iwp) :: i !< running index, x-dimension 2637 2731 INTEGER(iwp) :: j !< running index, y-dimension 2732 2638 2733 ! 2639 2734 !-- Dirichlet boundary conditions are required as the surface fluxes are … … 2680 2775 ENDIF 2681 2776 ! 2682 !-- Check if building types are set within a valid range. 2683 IF ( building_type < LBOUND( building_pars, 2 ) .AND. &2777 !-- Check if building types are set within a valid range. 2778 IF ( building_type < LBOUND( building_pars, 2 ) .AND. & 2684 2779 building_type > UBOUND( building_pars, 2 ) ) THEN 2685 2780 WRITE( message_string, * ) 'building_type = ', building_type, & … … 2696 2791 WRITE( message_string, * ) 'building_type = is out of ' // & 2697 2792 'the valid range at (j,i) = ', j, i 2698 CALL message( 'usm_check_parameters', 'PA0529', 2, 2, 0, 6, 0 ) 2793 CALL message( 'usm_check_parameters', 'PA0529', 2, 2, 0, 6, 0 ) 2699 2794 ENDIF 2700 2795 ENDDO 2701 2796 ENDDO 2702 2797 ENDIF 2703 2704 2798 END SUBROUTINE usm_check_parameters 2705 2799 … … 2719 2813 IMPLICIT NONE 2720 2814 2721 INTEGER(iwp), INTENT(IN) :: av !< 2722 CHARACTER (len=*), INTENT(IN) :: variable !< 2815 INTEGER(iwp), INTENT(IN) :: av !< flag if averaged 2816 CHARACTER (len=*), INTENT(IN) :: variable !< variable name 2723 2817 INTEGER(iwp), INTENT(IN) :: nzb_do !< lower limit of the data output (usually 0) 2724 2818 INTEGER(iwp), INTENT(IN) :: nzt_do !< vertical upper limit of the data output (usually nz_do3d) … … 2727 2821 REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) :: temp_pf !< temp array for urban surface output procedure 2728 2822 2729 CHARACTER (len=varnamelength) :: var 2730 INTEGER(iwp), PARAMETER :: nd = 5 2823 CHARACTER (len=varnamelength) :: var !< trimmed variable name 2824 INTEGER(iwp), PARAMETER :: nd = 5 !< number of directions 2731 2825 CHARACTER(len=6), DIMENSION(0:nd-1), PARAMETER :: dirname = (/ '_roof ', '_south', '_north', '_west ', '_east ' /) 2732 2826 INTEGER(iwp), DIMENSION(0:nd-1), PARAMETER :: dirint = (/ iup_u, isouth_u, inorth_u, iwest_u, ieast_u /) … … 2734 2828 !< index for surf_*_v: 0:3 = (North, South, East, West) 2735 2829 INTEGER(iwp) :: ids,idsint,idsidx,isvf 2736 INTEGER(iwp) :: i,j,k,iwl,istat, l, m 2830 INTEGER(iwp) :: i,j,k,iwl,istat, l, m !< running indices 2737 2831 2738 2832 found = .TRUE. … … 2756 2850 ENDIF 2757 2851 IF ( var(1:11) == 'usm_t_wall_' .AND. len(TRIM(var)) >= 12 ) THEN 2852 ! 2758 2853 !-- wall layers 2759 2854 READ(var(12:12), '(I1)', iostat=istat ) iwl … … 2763 2858 ENDIF 2764 2859 IF ( var(1:13) == 'usm_t_window_' .AND. len(TRIM(var)) >= 14 ) THEN 2860 ! 2765 2861 !-- window layers 2766 2862 READ(var(14:14), '(I1)', iostat=istat ) iwl … … 2770 2866 ENDIF 2771 2867 IF ( var(1:12) == 'usm_t_green_' .AND. len(TRIM(var)) >= 13 ) THEN 2868 ! 2772 2869 !-- green layers 2773 2870 READ(var(13:13), '(I1)', iostat=istat ) iwl … … 2777 2874 ENDIF 2778 2875 IF ( var(1:8) == 'usm_swc_' .AND. len(TRIM(var)) >= 9 ) THEN 2876 ! 2779 2877 !-- green layers soil water content 2780 2878 READ(var(9:9), '(I1)', iostat=istat ) iwl … … 2787 2885 2788 2886 CASE ( 'usm_surfz' ) 2887 ! 2789 2888 !-- array of surface height (z) 2790 2889 IF ( idsint == iup_u ) THEN … … 2806 2905 2807 2906 CASE ( 'usm_surfcat' ) 2907 ! 2808 2908 !-- surface category 2809 2909 IF ( idsint == iup_u ) THEN … … 2825 2925 2826 2926 CASE ( 'usm_surfalb' ) 2927 ! 2827 2928 !-- surface albedo, weighted average 2828 2929 IF ( idsint == iup_u ) THEN … … 2854 2955 2855 2956 CASE ( 'usm_surfemis' ) 2957 ! 2856 2958 !-- surface emissivity, weighted average 2857 2959 IF ( idsint == iup_u ) THEN … … 2883 2985 2884 2986 CASE ( 'usm_surfwintrans' ) 2987 ! 2885 2988 !-- transmissivity window tiles 2886 2989 IF ( idsint == iup_u ) THEN … … 2902 3005 2903 3006 CASE ( 'usm_wshf' ) 3007 ! 2904 3008 !-- array of sensible heat flux from surfaces 2905 3009 IF ( av == 0 ) THEN … … 2941 3045 2942 3046 CASE ( 'usm_qsws' ) 3047 ! 2943 3048 !-- array of latent heat flux from surfaces 2944 3049 IF ( av == 0 ) THEN … … 2979 3084 2980 3085 CASE ( 'usm_qsws_veg' ) 3086 ! 2981 3087 !-- array of latent heat flux from vegetation surfaces 2982 3088 IF ( av == 0 ) THEN … … 3017 3123 3018 3124 CASE ( 'usm_qsws_liq' ) 3125 ! 3019 3126 !-- array of latent heat flux from surfaces with liquid 3020 3127 IF ( av == 0 ) THEN … … 3055 3162 3056 3163 CASE ( 'usm_wghf' ) 3164 ! 3057 3165 !-- array of heat flux from ground (land, wall, roof) 3058 3166 IF ( av == 0 ) THEN … … 3093 3201 3094 3202 CASE ( 'usm_wghf_window' ) 3203 ! 3095 3204 !-- array of heat flux from window ground (land, wall, roof) 3096 3097 3205 IF ( av == 0 ) THEN 3098 3206 IF ( idsint == iup_u ) THEN … … 3132 3240 3133 3241 CASE ( 'usm_wghf_green' ) 3242 ! 3134 3243 !-- array of heat flux from green ground (land, wall, roof) 3135 3136 3244 IF ( av == 0 ) THEN 3137 3245 IF ( idsint == iup_u ) THEN … … 3171 3279 3172 3280 CASE ( 'usm_iwghf' ) 3281 ! 3173 3282 !-- array of heat flux from indoor ground (land, wall, roof) 3174 3283 IF ( av == 0 ) THEN … … 3209 3318 3210 3319 CASE ( 'usm_iwghf_window' ) 3320 ! 3211 3321 !-- array of heat flux from indoor window ground (land, wall, roof) 3212 3213 3322 IF ( av == 0 ) THEN 3214 3323 IF ( idsint == iup_u ) THEN … … 3248 3357 3249 3358 CASE ( 'usm_t_surf_wall' ) 3359 ! 3250 3360 !-- surface temperature for surfaces 3251 3361 IF ( av == 0 ) THEN … … 3286 3396 3287 3397 CASE ( 'usm_t_surf_window' ) 3398 ! 3288 3399 !-- surface temperature for window surfaces 3289 3290 3400 IF ( av == 0 ) THEN 3291 3401 IF ( idsint == iup_u ) THEN … … 3328 3438 3329 3439 CASE ( 'usm_t_surf_green' ) 3440 ! 3330 3441 !-- surface temperature for green surfaces 3331 3332 3442 IF ( av == 0 ) THEN 3333 3443 IF ( idsint == iup_u ) THEN … … 3370 3480 3371 3481 CASE ( 'usm_theta_10cm' ) 3482 ! 3372 3483 !-- near surface temperature for whole surfaces 3373 3374 3484 IF ( av == 0 ) THEN 3375 3485 IF ( idsint == iup_u ) THEN … … 3412 3522 3413 3523 CASE ( 'usm_t_wall' ) 3524 ! 3414 3525 !-- wall temperature for iwl layer of walls and land 3415 3526 IF ( av == 0 ) THEN … … 3450 3561 3451 3562 CASE ( 'usm_t_window' ) 3563 ! 3452 3564 !-- window temperature for iwl layer of walls and land 3453 3565 IF ( av == 0 ) THEN … … 3488 3600 3489 3601 CASE ( 'usm_t_green' ) 3602 ! 3490 3603 !-- green temperature for iwl layer of walls and land 3491 3604 IF ( av == 0 ) THEN … … 3526 3639 3527 3640 CASE ( 'usm_swc' ) 3641 ! 3528 3642 !-- soil water content for iwl layer of walls and land 3529 3643 IF ( av == 0 ) THEN … … 3612 3726 var(1:16) == 'usm_t_surf_green' .OR. var(1:11) == 'usm_t_green' .OR. & 3613 3727 var(1:15) == 'usm_theta_10cm' .OR. & 3614 var(1:9) == 'usm_surfz' .OR. var(1:11) == 'usm_surfcat' .OR. 3728 var(1:9) == 'usm_surfz' .OR. var(1:11) == 'usm_surfcat' .OR. & 3615 3729 var(1:11) == 'usm_surfalb' .OR. var(1:12) == 'usm_surfemis' .OR. & 3616 3730 var(1:16) == 'usm_surfwintrans' .OR. var(1:7) == 'usm_swc' ) THEN … … 3642 3756 3643 3757 CALL location_message( ' initialization of wall surface model', .TRUE. ) 3644 3758 3759 ! 3645 3760 !-- Calculate wall grid spacings. 3646 3761 !-- Temperature is defined at the center of the wall layers, … … 3659 3774 surf_usm_h%zw_window(k-1,m) 3660 3775 ENDDO 3661 ! surf_usm_h%dz_green(nzb_wall,m) = surf_usm_h%zw_green(nzb_wall,m)3662 ! DO k = nzb_wall+1, nzt_wall3663 ! surf_usm_h%dz_green(k,m) = surf_usm_h%zw_green(k,m) - &3664 ! surf_usm_h%zw_green(k-1,m)3665 ! ENDDO3666 3776 3667 3777 surf_usm_h%dz_wall(nzt_wall+1,m) = surf_usm_h%dz_wall(nzt_wall,m) … … 3681 3791 surf_usm_h%dz_window_stag(nzt_wall,m) = surf_usm_h%dz_window(nzt_wall,m) 3682 3792 3683 ! surf_usm_h%dz_green(nzt_wall+1,m) = surf_usm_h%dz_green(nzt_wall,m) 3684 ! 3685 ! DO k = nzb_wall, nzt_wall-1 3686 ! surf_usm_h%dz_green_stag(k,m) = 0.5 * ( & 3687 ! surf_usm_h%dz_green(k+1,m) + surf_usm_h%dz_green(k,m) ) 3688 ! ENDDO 3689 ! surf_usm_h%dz_green_stag(nzt_wall,m) = surf_usm_h%dz_green(nzt_wall,m) 3690 !------------- 3691 IF (surf_usm_h%green_type_roof(m) == 2.0_wp ) then 3692 soil_type = 3 !extensiv green roof 3693 surf_usm_h%lai(m) = 2.0_wp 3793 IF (surf_usm_h%green_type_roof(m) == 2.0_wp ) THEN 3794 ! 3795 !-- extensive green roof 3796 !-- set ratio of substrate layer thickness, soil-type and LAI 3797 soil_type = 3 3798 surf_usm_h%lai(m) = 2.0_wp 3694 3799 3695 surf_usm_h%zw_green(nzb_wall,m)= 0.05_wp3696 surf_usm_h%zw_green(nzb_wall+1,m)= 0.10_wp3697 surf_usm_h%zw_green(nzb_wall+2,m)= 0.15_wp3698 surf_usm_h%zw_green(nzb_wall+3,m)= 0.20_wp3800 surf_usm_h%zw_green(nzb_wall,m) = 0.05_wp 3801 surf_usm_h%zw_green(nzb_wall+1,m) = 0.10_wp 3802 surf_usm_h%zw_green(nzb_wall+2,m) = 0.15_wp 3803 surf_usm_h%zw_green(nzb_wall+3,m) = 0.20_wp 3699 3804 ELSE 3700 soil_type = 6 !intensiv green roof 3701 surf_usm_h%lai(m) = 4.0_wp 3805 ! 3806 !-- intensiv green roof 3807 !-- set ratio of substrate layer thickness, soil-type and LAI 3808 soil_type = 6 3809 surf_usm_h%lai(m) = 4.0_wp 3702 3810 3703 surf_usm_h%zw_green(nzb_wall,m)= 0.05_wp3704 surf_usm_h%zw_green(nzb_wall+1,m)= 0.10_wp3705 surf_usm_h%zw_green(nzb_wall+2,m)= 0.40_wp3706 surf_usm_h%zw_green(nzb_wall+3,m)= 0.80_wp3811 surf_usm_h%zw_green(nzb_wall,m) = 0.05_wp 3812 surf_usm_h%zw_green(nzb_wall+1,m) = 0.10_wp 3813 surf_usm_h%zw_green(nzb_wall+2,m) = 0.40_wp 3814 surf_usm_h%zw_green(nzb_wall+3,m) = 0.80_wp 3707 3815 ENDIF 3708 3816 … … 3758 3866 surf_usm_h%l_vg_green(m) = l_vangenuchten 3759 3867 surf_usm_h%n_vg_green(m) = n_vangenuchten 3760 surf_usm_h%gamma_w_green_sat(k,m) 3868 surf_usm_h%gamma_w_green_sat(k,m) = hydraulic_conductivity 3761 3869 swc_sat_h(k,m) = saturation_moisture 3762 3870 fc_h(k,m) = field_capacity … … 3764 3872 swc_res_h(k,m) = residual_moisture 3765 3873 ENDDO 3766 !------------------------------- 3874 3767 3875 ENDDO 3876 3768 3877 surf_usm_h%ddz_wall = 1.0_wp / surf_usm_h%dz_wall 3769 3878 surf_usm_h%ddz_wall_stag = 1.0_wp / surf_usm_h%dz_wall_stag … … 3802 3911 surf_usm_v(l)%dz_wall_stag(nzt_wall,m) = & 3803 3912 surf_usm_v(l)%dz_wall(nzt_wall,m) 3804 surf_usm_v(l)%dz_window(nzt_wall+1,m) = 3913 surf_usm_v(l)%dz_window(nzt_wall+1,m) = & 3805 3914 surf_usm_v(l)%dz_window(nzt_wall,m) 3806 3915 … … 3812 3921 surf_usm_v(l)%dz_window_stag(nzt_wall,m) = & 3813 3922 surf_usm_v(l)%dz_window(nzt_wall,m) 3814 surf_usm_v(l)%dz_green(nzt_wall+1,m) = &3923 surf_usm_v(l)%dz_green(nzt_wall+1,m) = & 3815 3924 surf_usm_v(l)%dz_green(nzt_wall,m) 3816 3925 … … 3831 3940 ENDDO 3832 3941 3833 3834 ! soil_type = 63835 ! !-- Initialize standard soil types. It is possible to overwrite each3836 ! !-- parameter by setting the respecticy NAMELIST variable to a3837 ! !-- value /= 9999999.9.3838 ! IF ( soil_type /= 0 ) THEN3839 !3840 ! IF ( alpha_vangenuchten == 9999999.9_wp ) THEN3841 ! alpha_vangenuchten = soil_pars(0,soil_type)3842 ! ENDIF3843 !3844 ! IF ( l_vangenuchten == 9999999.9_wp ) THEN3845 ! l_vangenuchten = soil_pars(1,soil_type)3846 ! ENDIF3847 !3848 ! IF ( n_vangenuchten == 9999999.9_wp ) THEN3849 ! n_vangenuchten = soil_pars(2,soil_type)3850 ! ENDIF3851 !3852 ! IF ( hydraulic_conductivity == 9999999.9_wp ) THEN3853 ! hydraulic_conductivity = soil_pars(3,soil_type)3854 ! ENDIF3855 !3856 ! IF ( saturation_moisture == 9999999.9_wp ) THEN3857 ! saturation_moisture = m_soil_pars(0,soil_type)3858 ! ENDIF3859 !3860 ! IF ( field_capacity == 9999999.9_wp ) THEN3861 ! field_capacity = m_soil_pars(1,soil_type)3862 ! ENDIF3863 !3864 ! IF ( wilting_point == 9999999.9_wp ) THEN3865 ! wilting_point = m_soil_pars(2,soil_type)3866 ! ENDIF3867 !3868 ! IF ( residual_moisture == 9999999.9_wp ) THEN3869 ! residual_moisture = m_soil_pars(3,soil_type)3870 ! ENDIF3871 !3872 ! DO m = 1, surf_usm_h%ns3873 ! DO k = nzb_wall, nzt_wall+13874 ! swc_h(k,m) = field_capacity3875 ! rootfr_h(k,m) = 0.5_wp3876 ! ENDDO3877 ! ENDDO3878 ! ! !3879 ! ! !-- Vertical surfaces3880 ! ! DO l = 0, 33881 ! ! DO m = 1, surf_usm_v(l)%ns3882 ! ! DO k = nzb_wall, nzt_wall+13883 ! ! swc_v(l)%t(k,m) = 0.5_wp3884 ! ! ENDDO3885 ! ! ENDDO3886 ! ! ENDDO3887 !3888 ! ENDIF3889 !3890 ! !3891 ! !-- Map values to the respective 2D arrays3892 ! surf_usm_h%alpha_vg_green = alpha_vangenuchten3893 ! surf_usm_h%l_vg_green = l_vangenuchten3894 ! surf_usm_h%n_vg_green = n_vangenuchten3895 ! surf_usm_h%gamma_w_green_sat = hydraulic_conductivity3896 ! swc_sat_h = saturation_moisture3897 ! fc_h = field_capacity3898 ! wilt_h = wilting_point3899 ! swc_res_h = residual_moisture3900 ! ! r_soil_min = min_soil_resistance3901 3942 3902 3943 CALL location_message( ' wall structures filed out', .TRUE. ) … … 3971 4012 3972 4013 CALL cpu_log( log_point_s(78), 'usm_init', 'start' ) 4014 ! 3973 4015 !-- surface forcing have to be disabled for LSF 3974 4016 !-- in case of enabled urban surface module … … 4037 4079 ! 4038 4080 !-- Set flag for ground level 4039 IF ( z_agl <= ground_floor_level_l ) 4081 IF ( z_agl <= ground_floor_level_l ) & 4040 4082 surf_usm_v(l)%ground_level(m) = .TRUE. 4041 4083 … … 4057 4099 ENDDO 4058 4100 4059 !---------------------------------------------------------------------------------------------4060 4101 ! 4061 4102 !-- Map values onto horizontal elemements 4062 4103 DO m = 1, surf_usm_h%ns 4063 surf_usm_h%r_canopy_min(m) = 200.0_wp ! min_canopy_resistance4064 surf_usm_h%g_d(m) = 0.0_wp !canopy_resistance_coefficient4104 surf_usm_h%r_canopy_min(m) = 200.0_wp !< min_canopy_resistance 4105 surf_usm_h%g_d(m) = 0.0_wp !< canopy_resistance_coefficient 4065 4106 ENDDO 4066 4107 ! … … 4069 4110 DO l = 0, 3 4070 4111 DO m = 1, surf_usm_v(l)%ns 4071 surf_usm_v(l)%r_canopy_min(m) = 200.0_wp ! min_canopy_resistance4072 surf_usm_v(l)%g_d(m) = 0.0_wp !canopy_resistance_coefficient4112 surf_usm_v(l)%r_canopy_min(m) = 200.0_wp !< min_canopy_resistance 4113 surf_usm_v(l)%g_d(m) = 0.0_wp !< canopy_resistance_coefficient 4073 4114 ENDDO 4074 4115 ENDDO 4075 !--------------------------------------------------------------------------------------------- 4076 ! 4116 4077 4117 ! 4078 4118 !-- Initialize urban-type surface attribute. According to initialization in … … 4637 4677 IF ( building_pars_f%pars_xy(ind_wall_frac_r,j,i) /= building_pars_f%fill ) & 4638 4678 surf_usm_h%frac(ind_veg_wall,m) = building_pars_f%pars_xy(ind_wall_frac_r,j,i) 4639 IF ( building_pars_f%pars_xy(ind_green_frac_r,j,i) /= building_pars_f%fill ) &4679 IF ( building_pars_f%pars_xy(ind_green_frac_r,j,i) /= building_pars_f%fill ) & 4640 4680 surf_usm_h%frac(ind_pav_green,m) = building_pars_f%pars_xy(ind_green_frac_r,j,i) 4641 4681 IF ( building_pars_f%pars_xy(ind_win_frac_r,j,i) /= building_pars_f%fill ) & … … 4643 4683 4644 4684 4645 IF ( building_pars_f%pars_xy(ind_lai_r,j,i) /= building_pars_f%fill ) &4685 IF ( building_pars_f%pars_xy(ind_lai_r,j,i) /= building_pars_f%fill ) & 4646 4686 surf_usm_h%lai(m) = building_pars_f%pars_xy(ind_lai_r,j,i) 4647 4687 … … 4720 4760 IF ( building_pars_f%pars_xy(ind_alb_wall_r,j,i) /= building_pars_f%fill ) & 4721 4761 surf_usm_h%albedo_type(ind_veg_wall,m) = building_pars_f%pars_xy(ind_alb_wall_r,j,i) 4722 IF ( building_pars_f%pars_xy(ind_alb_green_r,j,i) /= building_pars_f%fill ) 4762 IF ( building_pars_f%pars_xy(ind_alb_green_r,j,i) /= building_pars_f%fill ) & 4723 4763 surf_usm_h%albedo_type(ind_pav_green,m) = building_pars_f%pars_xy(ind_alb_green_r,j,i) 4724 IF ( building_pars_f%pars_xy(ind_alb_win_r,j,i) /= building_pars_f%fill ) &4764 IF ( building_pars_f%pars_xy(ind_alb_win_r,j,i) /= building_pars_f%fill ) & 4725 4765 surf_usm_h%albedo_type(ind_wat_win,m) = building_pars_f%pars_xy(ind_alb_win_r,j,i) 4726 4766 … … 4741 4781 IF ( building_pars_f%pars_xy(ind_thick_4_wall_r,j,i) /= building_pars_f%fill ) & 4742 4782 surf_usm_h%zw_green(nzb_wall+3,m) = building_pars_f%pars_xy(ind_thick_4_wall_r,j,i) 4743 IF ( building_pars_f%pars_xy(ind_thick_1_win_r,j,i) /= building_pars_f%fill ) &4783 IF ( building_pars_f%pars_xy(ind_thick_1_win_r,j,i) /= building_pars_f%fill ) & 4744 4784 surf_usm_h%zw_window(nzb_wall,m) = building_pars_f%pars_xy(ind_thick_1_win_r,j,i) 4745 IF ( building_pars_f%pars_xy(ind_thick_2_win_r,j,i) /= building_pars_f%fill ) &4785 IF ( building_pars_f%pars_xy(ind_thick_2_win_r,j,i) /= building_pars_f%fill ) & 4746 4786 surf_usm_h%zw_window(nzb_wall+1,m) = building_pars_f%pars_xy(ind_thick_2_win_r,j,i) 4747 IF ( building_pars_f%pars_xy(ind_thick_3_win_r,j,i) /= building_pars_f%fill ) &4787 IF ( building_pars_f%pars_xy(ind_thick_3_win_r,j,i) /= building_pars_f%fill ) & 4748 4788 surf_usm_h%zw_window(nzb_wall+2,m) = building_pars_f%pars_xy(ind_thick_3_win_r,j,i) 4749 IF ( building_pars_f%pars_xy(ind_thick_4_win_r,j,i) /= building_pars_f%fill ) &4789 IF ( building_pars_f%pars_xy(ind_thick_4_win_r,j,i) /= building_pars_f%fill ) & 4750 4790 surf_usm_h%zw_window(nzb_wall+3,m) = building_pars_f%pars_xy(ind_thick_4_win_r,j,i) 4751 4791 … … 4783 4823 ind_alb_win = MERGE( ind_alb_win_gfl, ind_alb_win_agfl, & 4784 4824 surf_usm_v(l)%ground_level(m) ) 4785 ind_wall_frac = MERGE( ind_wall_frac_gfl, ind_wall_frac_agfl, 4825 ind_wall_frac = MERGE( ind_wall_frac_gfl, ind_wall_frac_agfl, & 4786 4826 surf_usm_v(l)%ground_level(m) ) 4787 4827 ind_win_frac = MERGE( ind_win_frac_gfl, ind_win_frac_agfl, & … … 4947 4987 IF ( building_pars_f%pars_xy(ind_alb_wall,j,i) /= building_pars_f%fill ) & 4948 4988 surf_usm_v(l)%albedo_type(ind_veg_wall,m) = building_pars_f%pars_xy(ind_alb_wall,j,i) 4949 IF ( building_pars_f%pars_xy(ind_alb_green,j,i) /= building_pars_f%fill ) 4989 IF ( building_pars_f%pars_xy(ind_alb_green,j,i) /= building_pars_f%fill ) & 4950 4990 surf_usm_v(l)%albedo_type(ind_pav_green,m) = building_pars_f%pars_xy(ind_alb_green,j,i) 4951 IF ( building_pars_f%pars_xy(ind_alb_win,j,i) /= building_pars_f%fill ) &4991 IF ( building_pars_f%pars_xy(ind_alb_win,j,i) /= building_pars_f%fill ) & 4952 4992 surf_usm_v(l)%albedo_type(ind_wat_win,m) = building_pars_f%pars_xy(ind_alb_win,j,i) 4953 4993 … … 5001 5041 CALL usm_read_urban_surface_types() 5002 5042 5003 !-- init material heat model5004 5043 CALL usm_init_material_model() 5005 5044 ! 5006 5045 !-- init anthropogenic sources of heat 5007 5046 IF ( usm_anthropogenic_heat ) THEN 5047 ! 5008 5048 !-- init anthropogenic sources of heat (from transportation for now) 5009 5049 CALL usm_read_anthropogenic_heat() … … 5046 5086 surf_usm_v(l)%z0(m) = 0.9_wp * surf_usm_v(l)%z_mo(m) 5047 5087 5048 WRITE( message_string, * ) 'z0 exceeds surface-layer height '// &5049 'at vertical urban surface and is ' // &5050 'decreased appropriately at grid point (i,j) = ', &5051 surf_usm_v(l)%i(m)+surf_usm_v(l)%ioff, &5088 WRITE( message_string, * ) 'z0 exceeds surface-layer height '// & 5089 'at vertical urban surface and is ' // & 5090 'decreased appropriately at grid point (i,j) = ', & 5091 surf_usm_v(l)%i(m)+surf_usm_v(l)%ioff, & 5052 5092 surf_usm_v(l)%j(m)+surf_usm_v(l)%joff 5053 CALL message( 'urban_surface_model_mod', 'PA0503', &5093 CALL message( 'urban_surface_model_mod', 'PA0503', & 5054 5094 0, 0, 0, 6, 0 ) 5055 5095 ENDIF … … 5059 5099 surf_usm_v(l)%z0q(m) = 0.9_wp * surf_usm_v(l)%z_mo(m) 5060 5100 5061 WRITE( message_string, * ) 'z0h exceeds surface-layer height '// &5062 'at vertical urban surface and is ' // &5063 'decreased appropriately at grid point (i,j) = ', &5064 surf_usm_v(l)%i(m)+surf_usm_v(l)%ioff, &5101 WRITE( message_string, * ) 'z0h exceeds surface-layer height '// & 5102 'at vertical urban surface and is ' // & 5103 'decreased appropriately at grid point (i,j) = ', & 5104 surf_usm_v(l)%i(m)+surf_usm_v(l)%ioff, & 5065 5105 surf_usm_v(l)%j(m)+surf_usm_v(l)%joff 5066 CALL message( 'urban_surface_model_mod', 'PA0507', &5106 CALL message( 'urban_surface_model_mod', 'PA0507', & 5067 5107 0, 0, 0, 6, 0 ) 5068 5108 ENDIF 5069 5109 ENDDO 5070 ENDDO 5071 5110 ENDDO 5111 5112 ! 5072 5113 !-- Intitialization of the surface and wall/ground/roof temperature 5073 5114 ! 5074 5115 !-- Initialization for restart runs 5075 5116 IF ( TRIM( initializing_actions ) /= 'read_restart_data' .AND. & … … 5117 5158 ENDDO 5118 5159 ENDIF 5119 5160 ! 5120 5161 !-- initial values for t_wall 5121 5162 !-- outer value is set to surface temperature … … 5172 5213 ENDIF 5173 5214 5215 ! 5174 5216 !-- If specified, replace constant wall temperatures with fully 3D values from file 5175 5217 IF ( read_wall_temp_3d ) CALL usm_read_wall_temperature() … … 5179 5221 CALL user_init_urban_surface 5180 5222 5223 ! 5181 5224 !-- initialize prognostic values for the first timestep 5182 5225 t_surf_wall_h_p = t_surf_wall_h … … 5194 5237 t_green_v_p = t_green_v 5195 5238 5239 ! 5196 5240 !-- Adjust radiative fluxes for urban surface at model start 5197 5241 !CALL radiation_interaction … … 5201 5245 m_liq_usm_h_p = m_liq_usm_h 5202 5246 m_liq_usm_v_p = m_liq_usm_v 5203 !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++5204 5247 ! 5205 5248 !-- Set initial values for prognostic quantities … … 5230 5273 ENDDO 5231 5274 ENDIF 5232 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 5233 5234 CALL cpu_log( log_point_s(78), 'usm_init', 'stop' ) 5235 5236 CALL location_message( 'finished', .TRUE. ) 5275 5276 CALL cpu_log( log_point_s(78), 'usm_init', 'stop' ) 5277 5278 CALL location_message( 'finished', .TRUE. ) 5237 5279 5238 5280 END SUBROUTINE usm_init … … 5243 5285 ! ------------ 5244 5286 ! 5245 !> Wall model as part of the urban surface model. The model predicts wall 5246 !> temperature. 5287 !> Wall model as part of the urban surface model. The model predicts vertical 5288 !> and horizontal wall / roof temperatures and window layer temperatures. 5289 !> No window layer temperature calculactions during spinup to increase 5290 !> possible timestep. 5247 5291 !------------------------------------------------------------------------------! 5248 5292 SUBROUTINE usm_material_heat_model( spinup ) … … 5254 5298 5255 5299 REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: wtend, wintend !< tendency 5256 REAL(wp) :: win_absorp !absorption coefficient from transmissivity5300 REAL(wp) :: win_absorp !< absorption coefficient from transmissivity 5257 5301 REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: wall_mod 5258 5302 5259 LOGICAL :: spinup !if true, no calculation of window temperatures5303 LOGICAL :: spinup !< if true, no calculation of window temperatures 5260 5304 5261 5305 wall_mod=1.0_wp 5262 if (usm_wall_mod .AND. spinup) then5263 dokw=nzb_wall,nzb_wall+15264 wall_mod(kw)=0.1_wp5265 enddo5266 endif5306 IF (usm_wall_mod .AND. spinup) THEN 5307 DO kw=nzb_wall,nzb_wall+1 5308 wall_mod(kw)=0.1_wp 5309 ENDDO 5310 ENDIF 5267 5311 5268 5312 ! … … 5291 5335 * ( surf_usm_h%lambda_h_green(nzt_wall,m)* wall_mod(nzt_wall) & 5292 5336 * surf_usm_h%ddz_green(nzt_wall,m) & 5293 + surf_usm_h%lambda_h(nzb_wall,m) * wall_mod(nzb_wall) 5337 + surf_usm_h%lambda_h(nzb_wall,m) * wall_mod(nzb_wall) & 5294 5338 * surf_usm_h%ddz_wall(nzb_wall,m) ) & 5295 5339 / ( surf_usm_h%ddz_green(nzt_wall,m) & … … 5298 5342 - t_green_h(nzt_wall,m) ) ) * & 5299 5343 surf_usm_h%ddz_wall_stag(nzb_wall,m) 5300 5301 IF ( indoor_model ) then 5344 ! 5345 !-- if indoor model ist used inner wall layer is calculated by using iwghf (indoor wall ground heat flux) 5346 IF ( indoor_model ) THEN 5302 5347 DO kw = nzb_wall+1, nzt_wall-1 5303 5348 wtend(kw) = (1.0_wp / surf_usm_h%rho_c_wall(kw,m)) & 5304 * ( surf_usm_h%lambda_h(kw,m) * wall_mod(kw) 5349 * ( surf_usm_h%lambda_h(kw,m) * wall_mod(kw) & 5305 5350 * ( t_wall_h(kw+1,m) - t_wall_h(kw,m) ) & 5306 5351 * surf_usm_h%ddz_wall(kw+1,m) & 5307 - surf_usm_h%lambda_h(kw-1,m) * wall_mod(kw-1) 5352 - surf_usm_h%lambda_h(kw-1,m) * wall_mod(kw-1) & 5308 5353 * ( t_wall_h(kw,m) - t_wall_h(kw-1,m) ) & 5309 5354 * surf_usm_h%ddz_wall(kw,m) & … … 5311 5356 ENDDO 5312 5357 wtend(nzt_wall) = (1.0_wp / surf_usm_h%rho_c_wall(nzt_wall,m)) * & 5313 ( -surf_usm_h%lambda_h(nzt_wall-1,m) * wall_mod(nzt_wall-1) * 5358 ( -surf_usm_h%lambda_h(nzt_wall-1,m) * wall_mod(nzt_wall-1) * & 5314 5359 ( t_wall_h(nzt_wall,m) & 5315 5360 - t_wall_h(nzt_wall-1,m) ) * & … … 5320 5365 DO kw = nzb_wall+1, nzt_wall 5321 5366 wtend(kw) = (1.0_wp / surf_usm_h%rho_c_wall(kw,m)) & 5322 * ( surf_usm_h%lambda_h(kw,m) * wall_mod(kw) 5367 * ( surf_usm_h%lambda_h(kw,m) * wall_mod(kw) & 5323 5368 * ( t_wall_h(kw+1,m) - t_wall_h(kw,m) ) & 5324 5369 * surf_usm_h%ddz_wall(kw+1,m) & 5325 - surf_usm_h%lambda_h(kw-1,m) * wall_mod(kw-1) 5370 - surf_usm_h%lambda_h(kw-1,m) * wall_mod(kw-1) & 5326 5371 * ( t_wall_h(kw,m) - t_wall_h(kw-1,m) ) & 5327 5372 * surf_usm_h%ddz_wall(kw,m) & … … 5335 5380 * surf_usm_h%tt_wall_m(nzb_wall:nzt_wall,m) ) 5336 5381 5337 if (.NOT. spinup) then 5338 win_absorp = -log(surf_usm_h%transmissivity(m)) / surf_usm_h%zw_window(nzt_wall,m) 5339 !-- prognostic equation for ground/roof window temperature t_window_h 5340 wintend(:) = 0.0_wp 5341 wintend(nzb_wall) = (1.0_wp / surf_usm_h%rho_c_window(nzb_wall,m)) * & 5342 ( surf_usm_h%lambda_h_window(nzb_wall,m) * & 5343 ( t_window_h(nzb_wall+1,m) & 5344 - t_window_h(nzb_wall,m) ) * & 5345 surf_usm_h%ddz_window(nzb_wall+1,m) & 5346 + surf_usm_h%wghf_eb_window(m) & 5347 + surf_usm_h%rad_sw_in(m) & 5348 * (1.0_wp - exp(-win_absorp & 5349 * surf_usm_h%zw_window(nzb_wall,m) ) ) & 5350 ) * surf_usm_h%ddz_window_stag(nzb_wall,m) 5351 5352 IF ( indoor_model ) then 5353 DO kw = nzb_wall+1, nzt_wall-1 5354 wintend(kw) = (1.0_wp / surf_usm_h%rho_c_window(kw,m)) & 5355 * ( surf_usm_h%lambda_h_window(kw,m) & 5356 * ( t_window_h(kw+1,m) - t_window_h(kw,m) ) & 5357 * surf_usm_h%ddz_window(kw+1,m) & 5358 - surf_usm_h%lambda_h_window(kw-1,m) & 5359 * ( t_window_h(kw,m) - t_window_h(kw-1,m) ) & 5360 * surf_usm_h%ddz_window(kw,m) & 5361 + surf_usm_h%rad_sw_in(m) & 5362 * (exp(-win_absorp & 5363 * surf_usm_h%zw_window(kw-1,m) ) & 5364 - exp(-win_absorp & 5365 * surf_usm_h%zw_window(kw,m) ) ) & 5366 ) * surf_usm_h%ddz_window_stag(kw,m) 5367 5368 ENDDO 5369 wintend(nzt_wall) = (1.0_wp / surf_usm_h%rho_c_window(nzt_wall,m)) * & 5370 ( -surf_usm_h%lambda_h_window(nzt_wall-1,m) * & 5371 ( t_window_h(nzt_wall,m) & 5372 - t_window_h(nzt_wall-1,m) ) * & 5373 surf_usm_h%ddz_window(nzt_wall,m) & 5374 + surf_usm_h%iwghf_eb_window(m) & 5375 + surf_usm_h%rad_sw_in(m) & 5376 * (exp(-win_absorp & 5377 * surf_usm_h%zw_window(nzt_wall-1,m) ) & 5378 - exp(-win_absorp & 5379 * surf_usm_h%zw_window(nzt_wall,m) ) ) & 5380 ) * surf_usm_h%ddz_window_stag(nzt_wall,m) 5381 ELSE 5382 DO kw = nzb_wall+1, nzt_wall 5383 wintend(kw) = (1.0_wp / surf_usm_h%rho_c_window(kw,m)) & 5384 * ( surf_usm_h%lambda_h_window(kw,m) & 5385 * ( t_window_h(kw+1,m) - t_window_h(kw,m) ) & 5386 * surf_usm_h%ddz_window(kw+1,m) & 5387 - surf_usm_h%lambda_h_window(kw-1,m) & 5388 * ( t_window_h(kw,m) - t_window_h(kw-1,m) ) & 5389 * surf_usm_h%ddz_window(kw,m) & 5390 + surf_usm_h%rad_sw_in(m) & 5391 * (exp(-win_absorp & 5392 * surf_usm_h%zw_window(kw-1,m) ) & 5393 - exp(-win_absorp & 5394 * surf_usm_h%zw_window(kw,m) ) ) & 5395 ) * surf_usm_h%ddz_window_stag(kw,m) 5396 5397 ENDDO 5398 ENDIF 5399 5400 t_window_h_p(nzb_wall:nzt_wall,m) = t_window_h(nzb_wall:nzt_wall,m) & 5382 ! 5383 !-- during spinup the tempeature inside window layers is not calculated to make larger timesteps possible 5384 IF ( .NOT. spinup) THEN 5385 win_absorp = -log(surf_usm_h%transmissivity(m)) / surf_usm_h%zw_window(nzt_wall,m) 5386 ! 5387 !-- prognostic equation for ground/roof window temperature t_window_h 5388 !-- takes absorption of shortwave radiation into account 5389 wintend(:) = 0.0_wp 5390 wintend(nzb_wall) = (1.0_wp / surf_usm_h%rho_c_window(nzb_wall,m)) * & 5391 ( surf_usm_h%lambda_h_window(nzb_wall,m) * & 5392 ( t_window_h(nzb_wall+1,m) & 5393 - t_window_h(nzb_wall,m) ) * & 5394 surf_usm_h%ddz_window(nzb_wall+1,m) & 5395 + surf_usm_h%wghf_eb_window(m) & 5396 + surf_usm_h%rad_sw_in(m) & 5397 * (1.0_wp - exp(-win_absorp & 5398 * surf_usm_h%zw_window(nzb_wall,m) ) ) & 5399 ) * surf_usm_h%ddz_window_stag(nzb_wall,m) 5400 5401 IF ( indoor_model ) THEN 5402 DO kw = nzb_wall+1, nzt_wall-1 5403 wintend(kw) = (1.0_wp / surf_usm_h%rho_c_window(kw,m)) & 5404 * ( surf_usm_h%lambda_h_window(kw,m) & 5405 * ( t_window_h(kw+1,m) - t_window_h(kw,m) ) & 5406 * surf_usm_h%ddz_window(kw+1,m) & 5407 - surf_usm_h%lambda_h_window(kw-1,m) & 5408 * ( t_window_h(kw,m) - t_window_h(kw-1,m) ) & 5409 * surf_usm_h%ddz_window(kw,m) & 5410 + surf_usm_h%rad_sw_in(m) & 5411 * (exp(-win_absorp & 5412 * surf_usm_h%zw_window(kw-1,m) ) & 5413 - exp(-win_absorp & 5414 * surf_usm_h%zw_window(kw,m) ) ) & 5415 ) * surf_usm_h%ddz_window_stag(kw,m) 5416 5417 ENDDO 5418 wintend(nzt_wall) = (1.0_wp / surf_usm_h%rho_c_window(nzt_wall,m)) * & 5419 ( -surf_usm_h%lambda_h_window(nzt_wall-1,m) * & 5420 ( t_window_h(nzt_wall,m) & 5421 - t_window_h(nzt_wall-1,m) ) * & 5422 surf_usm_h%ddz_window(nzt_wall,m) & 5423 + surf_usm_h%iwghf_eb_window(m) & 5424 + surf_usm_h%rad_sw_in(m) & 5425 * (exp(-win_absorp & 5426 * surf_usm_h%zw_window(nzt_wall-1,m) ) & 5427 - exp(-win_absorp & 5428 * surf_usm_h%zw_window(nzt_wall,m) ) ) & 5429 ) * surf_usm_h%ddz_window_stag(nzt_wall,m) 5430 ELSE 5431 DO kw = nzb_wall+1, nzt_wall 5432 wintend(kw) = (1.0_wp / surf_usm_h%rho_c_window(kw,m)) & 5433 * ( surf_usm_h%lambda_h_window(kw,m) & 5434 * ( t_window_h(kw+1,m) - t_window_h(kw,m) ) & 5435 * surf_usm_h%ddz_window(kw+1,m) & 5436 - surf_usm_h%lambda_h_window(kw-1,m) & 5437 * ( t_window_h(kw,m) - t_window_h(kw-1,m) ) & 5438 * surf_usm_h%ddz_window(kw,m) & 5439 + surf_usm_h%rad_sw_in(m) & 5440 * (exp(-win_absorp & 5441 * surf_usm_h%zw_window(kw-1,m) ) & 5442 - exp(-win_absorp & 5443 * surf_usm_h%zw_window(kw,m) ) ) & 5444 ) * surf_usm_h%ddz_window_stag(kw,m) 5445 5446 ENDDO 5447 ENDIF 5448 5449 t_window_h_p(nzb_wall:nzt_wall,m) = t_window_h(nzb_wall:nzt_wall,m) & 5401 5450 + dt_3d * ( tsc(2) & 5402 5451 * wintend(nzb_wall:nzt_wall) + tsc(3) & 5403 5452 * surf_usm_h%tt_window_m(nzb_wall:nzt_wall,m) ) 5404 5453 5405 endif 5454 ENDIF 5406 5455 5407 5456 ! … … 5421 5470 ENDIF 5422 5471 5423 if (.NOT. spinup) then 5424 !-- calculate t_window tendencies for the next Runge-Kutta step 5425 IF ( timestep_scheme(1:5) == 'runge' ) THEN 5426 IF ( intermediate_timestep_count == 1 ) THEN 5427 DO kw = nzb_wall, nzt_wall 5428 surf_usm_h%tt_window_m(kw,m) = wintend(kw) 5429 ENDDO 5430 ELSEIF ( intermediate_timestep_count < & 5431 intermediate_timestep_count_max ) THEN 5432 DO kw = nzb_wall, nzt_wall 5433 surf_usm_h%tt_window_m(kw,m) = -9.5625_wp * wintend(kw) + & 5434 5.3125_wp * surf_usm_h%tt_window_m(kw,m) 5435 ENDDO 5436 ENDIF 5472 IF (.NOT. spinup) THEN 5473 ! 5474 !-- calculate t_window tendencies for the next Runge-Kutta step 5475 IF ( timestep_scheme(1:5) == 'runge' ) THEN 5476 IF ( intermediate_timestep_count == 1 ) THEN 5477 DO kw = nzb_wall, nzt_wall 5478 surf_usm_h%tt_window_m(kw,m) = wintend(kw) 5479 ENDDO 5480 ELSEIF ( intermediate_timestep_count < & 5481 intermediate_timestep_count_max ) THEN 5482 DO kw = nzb_wall, nzt_wall 5483 surf_usm_h%tt_window_m(kw,m) = -9.5625_wp * wintend(kw) + & 5484 5.3125_wp * surf_usm_h%tt_window_m(kw,m) 5485 ENDDO 5486 ENDIF 5487 ENDIF 5437 5488 ENDIF 5438 5439 endif5440 5489 5441 5490 ENDDO … … 5454 5503 wtend(:) = 0.0_wp 5455 5504 5456 5457 5458 5459 5460 5461 5462 5463 5464 5465 5466 5467 5468 5469 5470 5471 5472 5473 5474 5475 5476 5477 5478 IF ( indoor_model ) then5505 wtend(nzb_wall) = (1.0_wp / surf_usm_v(l)%rho_c_wall(nzb_wall,m)) * & 5506 ( surf_usm_v(l)%lambda_h(nzb_wall,m) * wall_mod(nzb_wall) * & 5507 ( t_wall_v(l)%t(nzb_wall+1,m) & 5508 - t_wall_v(l)%t(nzb_wall,m) ) * & 5509 surf_usm_v(l)%ddz_wall(nzb_wall+1,m) & 5510 + surf_usm_v(l)%frac(ind_veg_wall,m) & 5511 / (surf_usm_v(l)%frac(ind_veg_wall,m) & 5512 + surf_usm_v(l)%frac(ind_pav_green,m) ) & 5513 * surf_usm_v(l)%wghf_eb(m) & 5514 - surf_usm_v(l)%frac(ind_pav_green,m) & 5515 / (surf_usm_v(l)%frac(ind_veg_wall,m) & 5516 + surf_usm_v(l)%frac(ind_pav_green,m) ) & 5517 * ( surf_usm_v(l)%lambda_h_green(nzt_wall,m)* wall_mod(nzt_wall) & 5518 * surf_usm_v(l)%ddz_green(nzt_wall,m) & 5519 + surf_usm_v(l)%lambda_h(nzb_wall,m)* wall_mod(nzb_wall) & 5520 * surf_usm_v(l)%ddz_wall(nzb_wall,m) ) & 5521 / ( surf_usm_v(l)%ddz_green(nzt_wall,m) & 5522 + surf_usm_v(l)%ddz_wall(nzb_wall,m) ) & 5523 * ( t_wall_v(l)%t(nzb_wall,m) & 5524 - t_green_v(l)%t(nzt_wall,m) ) ) * & 5525 surf_usm_v(l)%ddz_wall_stag(nzb_wall,m) 5526 5527 IF ( indoor_model ) THEN 5479 5528 DO kw = nzb_wall+1, nzt_wall-1 5480 5529 wtend(kw) = (1.0_wp / surf_usm_v(l)%rho_c_wall(kw,m)) & 5481 * ( surf_usm_v(l)%lambda_h(kw,m) * wall_mod(kw) 5530 * ( surf_usm_v(l)%lambda_h(kw,m) * wall_mod(kw) & 5482 5531 * ( t_wall_v(l)%t(kw+1,m) - t_wall_v(l)%t(kw,m) )& 5483 5532 * surf_usm_v(l)%ddz_wall(kw+1,m) & 5484 - surf_usm_v(l)%lambda_h(kw-1,m) * wall_mod(kw-1) 5533 - surf_usm_v(l)%lambda_h(kw-1,m) * wall_mod(kw-1) & 5485 5534 * ( t_wall_v(l)%t(kw,m) - t_wall_v(l)%t(kw-1,m) )& 5486 5535 * surf_usm_v(l)%ddz_wall(kw,m) & … … 5497 5546 DO kw = nzb_wall+1, nzt_wall 5498 5547 wtend(kw) = (1.0_wp / surf_usm_v(l)%rho_c_wall(kw,m)) & 5499 * ( surf_usm_v(l)%lambda_h(kw,m) * wall_mod(kw) 5548 * ( surf_usm_v(l)%lambda_h(kw,m) * wall_mod(kw) & 5500 5549 * ( t_wall_v(l)%t(kw+1,m) - t_wall_v(l)%t(kw,m) )& 5501 5550 * surf_usm_v(l)%ddz_wall(kw+1,m) & 5502 - surf_usm_v(l)%lambda_h(kw-1,m) * wall_mod(kw-1) 5551 - surf_usm_v(l)%lambda_h(kw-1,m) * wall_mod(kw-1) & 5503 5552 * ( t_wall_v(l)%t(kw,m) - t_wall_v(l)%t(kw-1,m) )& 5504 5553 * surf_usm_v(l)%ddz_wall(kw,m) & … … 5513 5562 * surf_usm_v(l)%tt_wall_m(nzb_wall:nzt_wall,m) ) 5514 5563 5515 if (.NOT. spinup) then 5516 win_absorp = -log(surf_usm_v(l)%transmissivity(m)) / surf_usm_v(l)%zw_window(nzt_wall,m) 5517 !-- prognostic equation for window temperature t_window_v 5518 wintend(:) = 0.0_wp 5519 wintend(nzb_wall) = (1.0_wp / surf_usm_v(l)%rho_c_window(nzb_wall,m)) * & 5520 ( surf_usm_v(l)%lambda_h_window(nzb_wall,m) * & 5521 ( t_window_v(l)%t(nzb_wall+1,m) & 5522 - t_window_v(l)%t(nzb_wall,m) ) * & 5523 surf_usm_v(l)%ddz_window(nzb_wall+1,m) & 5524 + surf_usm_v(l)%wghf_eb_window(m) & 5525 + surf_usm_v(l)%rad_sw_in(m) & 5526 * (1.0_wp - exp(-win_absorp & 5527 * surf_usm_v(l)%zw_window(nzb_wall,m) ) ) & 5528 ) * surf_usm_v(l)%ddz_window_stag(nzb_wall,m) 5529 5530 IF ( indoor_model ) then 5531 DO kw = nzb_wall+1, nzt_wall -1 5532 wintend(kw) = (1.0_wp / surf_usm_v(l)%rho_c_window(kw,m)) & 5533 * ( surf_usm_v(l)%lambda_h_window(kw,m) & 5534 * ( t_window_v(l)%t(kw+1,m) - t_window_v(l)%t(kw,m) ) & 5535 * surf_usm_v(l)%ddz_window(kw+1,m) & 5536 - surf_usm_v(l)%lambda_h_window(kw-1,m) & 5537 * ( t_window_v(l)%t(kw,m) - t_window_v(l)%t(kw-1,m) ) & 5538 * surf_usm_v(l)%ddz_window(kw,m) & 5539 + surf_usm_v(l)%rad_sw_in(m) & 5540 * (exp(-win_absorp & 5541 * surf_usm_v(l)%zw_window(kw-1,m) ) & 5542 - exp(-win_absorp & 5543 * surf_usm_v(l)%zw_window(kw,m) ) ) & 5544 ) * surf_usm_v(l)%ddz_window_stag(kw,m) 5545 ENDDO 5546 wintend(nzt_wall) = (1.0_wp / surf_usm_v(l)%rho_c_window(nzt_wall,m)) * & 5547 ( -surf_usm_v(l)%lambda_h_window(nzt_wall-1,m) * & 5548 ( t_window_v(l)%t(nzt_wall,m) & 5549 - t_window_v(l)%t(nzt_wall-1,m) ) * & 5550 surf_usm_v(l)%ddz_window(nzt_wall,m) & 5551 + surf_usm_v(l)%iwghf_eb_window(m) & 5552 + surf_usm_v(l)%rad_sw_in(m) & 5553 * (exp(-win_absorp & 5554 * surf_usm_v(l)%zw_window(nzt_wall-1,m) ) & 5555 - exp(-win_absorp & 5556 * surf_usm_v(l)%zw_window(nzt_wall,m) ) ) & 5557 ) * surf_usm_v(l)%ddz_window_stag(nzt_wall,m) 5558 ELSE 5559 DO kw = nzb_wall+1, nzt_wall 5560 wintend(kw) = (1.0_wp / surf_usm_v(l)%rho_c_window(kw,m)) & 5561 * ( surf_usm_v(l)%lambda_h_window(kw,m) & 5562 * ( t_window_v(l)%t(kw+1,m) - t_window_v(l)%t(kw,m) ) & 5563 * surf_usm_v(l)%ddz_window(kw+1,m) & 5564 - surf_usm_v(l)%lambda_h_window(kw-1,m) & 5565 * ( t_window_v(l)%t(kw,m) - t_window_v(l)%t(kw-1,m) ) & 5566 * surf_usm_v(l)%ddz_window(kw,m) & 5567 + surf_usm_v(l)%rad_sw_in(m) & 5568 * (exp(-win_absorp & 5569 * surf_usm_v(l)%zw_window(kw-1,m) ) & 5570 - exp(-win_absorp & 5571 * surf_usm_v(l)%zw_window(kw,m) ) ) & 5572 ) * surf_usm_v(l)%ddz_window_stag(kw,m) 5573 ENDDO 5564 IF (.NOT. spinup) THEN 5565 win_absorp = -log(surf_usm_v(l)%transmissivity(m)) / surf_usm_v(l)%zw_window(nzt_wall,m) 5566 ! 5567 !-- prognostic equation for window temperature t_window_v 5568 wintend(:) = 0.0_wp 5569 wintend(nzb_wall) = (1.0_wp / surf_usm_v(l)%rho_c_window(nzb_wall,m)) * & 5570 ( surf_usm_v(l)%lambda_h_window(nzb_wall,m) * & 5571 ( t_window_v(l)%t(nzb_wall+1,m) & 5572 - t_window_v(l)%t(nzb_wall,m) ) * & 5573 surf_usm_v(l)%ddz_window(nzb_wall+1,m) & 5574 + surf_usm_v(l)%wghf_eb_window(m) & 5575 + surf_usm_v(l)%rad_sw_in(m) & 5576 * (1.0_wp - exp(-win_absorp & 5577 * surf_usm_v(l)%zw_window(nzb_wall,m) ) ) & 5578 ) * surf_usm_v(l)%ddz_window_stag(nzb_wall,m) 5579 5580 IF ( indoor_model ) THEN 5581 DO kw = nzb_wall+1, nzt_wall -1 5582 wintend(kw) = (1.0_wp / surf_usm_v(l)%rho_c_window(kw,m)) & 5583 * ( surf_usm_v(l)%lambda_h_window(kw,m) & 5584 * ( t_window_v(l)%t(kw+1,m) - t_window_v(l)%t(kw,m) ) & 5585 * surf_usm_v(l)%ddz_window(kw+1,m) & 5586 - surf_usm_v(l)%lambda_h_window(kw-1,m) & 5587 * ( t_window_v(l)%t(kw,m) - t_window_v(l)%t(kw-1,m) ) & 5588 * surf_usm_v(l)%ddz_window(kw,m) & 5589 + surf_usm_v(l)%rad_sw_in(m) & 5590 * (exp(-win_absorp & 5591 * surf_usm_v(l)%zw_window(kw-1,m) ) & 5592 - exp(-win_absorp & 5593 * surf_usm_v(l)%zw_window(kw,m) ) ) & 5594 ) * surf_usm_v(l)%ddz_window_stag(kw,m) 5595 ENDDO 5596 wintend(nzt_wall) = (1.0_wp / surf_usm_v(l)%rho_c_window(nzt_wall,m)) * & 5597 ( -surf_usm_v(l)%lambda_h_window(nzt_wall-1,m) * & 5598 ( t_window_v(l)%t(nzt_wall,m) & 5599 - t_window_v(l)%t(nzt_wall-1,m) ) * & 5600 surf_usm_v(l)%ddz_window(nzt_wall,m) & 5601 + surf_usm_v(l)%iwghf_eb_window(m) & 5602 + surf_usm_v(l)%rad_sw_in(m) & 5603 * (exp(-win_absorp & 5604 * surf_usm_v(l)%zw_window(nzt_wall-1,m) ) & 5605 - exp(-win_absorp & 5606 * surf_usm_v(l)%zw_window(nzt_wall,m) ) ) & 5607 ) * surf_usm_v(l)%ddz_window_stag(nzt_wall,m) 5608 ELSE 5609 DO kw = nzb_wall+1, nzt_wall 5610 wintend(kw) = (1.0_wp / surf_usm_v(l)%rho_c_window(kw,m)) & 5611 * ( surf_usm_v(l)%lambda_h_window(kw,m) & 5612 * ( t_window_v(l)%t(kw+1,m) - t_window_v(l)%t(kw,m) ) & 5613 * surf_usm_v(l)%ddz_window(kw+1,m) & 5614 - surf_usm_v(l)%lambda_h_window(kw-1,m) & 5615 * ( t_window_v(l)%t(kw,m) - t_window_v(l)%t(kw-1,m) ) & 5616 * surf_usm_v(l)%ddz_window(kw,m) & 5617 + surf_usm_v(l)%rad_sw_in(m) & 5618 * (exp(-win_absorp & 5619 * surf_usm_v(l)%zw_window(kw-1,m) ) & 5620 - exp(-win_absorp & 5621 * surf_usm_v(l)%zw_window(kw,m) ) ) & 5622 ) * surf_usm_v(l)%ddz_window_stag(kw,m) 5623 ENDDO 5624 ENDIF 5625 5626 t_window_v_p(l)%t(nzb_wall:nzt_wall,m) = & 5627 t_window_v(l)%t(nzb_wall:nzt_wall,m) & 5628 + dt_3d * ( tsc(2) & 5629 * wintend(nzb_wall:nzt_wall) + tsc(3) & 5630 * surf_usm_v(l)%tt_window_m(nzb_wall:nzt_wall,m) ) 5574 5631 ENDIF 5575 5576 t_window_v_p(l)%t(nzb_wall:nzt_wall,m) = &5577 t_window_v(l)%t(nzb_wall:nzt_wall,m) &5578 + dt_3d * ( tsc(2) &5579 * wintend(nzb_wall:nzt_wall) + tsc(3) &5580 * surf_usm_v(l)%tt_window_m(nzb_wall:nzt_wall,m) )5581 endif5582 5632 5583 5633 ! … … 5599 5649 5600 5650 5601 if (.NOT. spinup) then 5602 !-- calculate t_window tendencies for the next Runge-Kutta step 5603 IF ( timestep_scheme(1:5) == 'runge' ) THEN 5604 IF ( intermediate_timestep_count == 1 ) THEN 5605 DO kw = nzb_wall, nzt_wall 5606 surf_usm_v(l)%tt_window_m(kw,m) = wintend(kw) 5607 ENDDO 5608 ELSEIF ( intermediate_timestep_count < & 5609 intermediate_timestep_count_max ) THEN 5610 DO kw = nzb_wall, nzt_wall 5611 surf_usm_v(l)%tt_window_m(kw,m) = & 5612 - 9.5625_wp * wintend(kw) + & 5613 5.3125_wp * surf_usm_v(l)%tt_window_m(kw,m) 5614 ENDDO 5615 ENDIF 5651 IF (.NOT. spinup) THEN 5652 ! 5653 !-- calculate t_window tendencies for the next Runge-Kutta step 5654 IF ( timestep_scheme(1:5) == 'runge' ) THEN 5655 IF ( intermediate_timestep_count == 1 ) THEN 5656 DO kw = nzb_wall, nzt_wall 5657 surf_usm_v(l)%tt_window_m(kw,m) = wintend(kw) 5658 ENDDO 5659 ELSEIF ( intermediate_timestep_count < & 5660 intermediate_timestep_count_max ) THEN 5661 DO kw = nzb_wall, nzt_wall 5662 surf_usm_v(l)%tt_window_m(kw,m) = & 5663 - 9.5625_wp * wintend(kw) + & 5664 5.3125_wp * surf_usm_v(l)%tt_window_m(kw,m) 5665 ENDDO 5666 ENDIF 5667 ENDIF 5616 5668 ENDIF 5617 endif5618 5669 5619 5670 ENDDO … … 5634 5685 IMPLICIT NONE 5635 5686 5636 INTEGER(iwp) :: i,j,k,l,kw, m 5637 5638 REAL(wp) :: ke, lambda_h_green_sat 5639 REAL(wp) :: h_vg !< Van Genuchten coef. h5640 REAL(wp) :: drho_l_lv 5687 INTEGER(iwp) :: i,j,k,l,kw, m !< running indices 5688 5689 REAL(wp) :: ke, lambda_h_green_sat !< heat conductivity for saturated soil 5690 REAL(wp) :: h_vg !< Van Genuchten coef. h 5691 REAL(wp) :: drho_l_lv !< frequently used parameter 5641 5692 5642 5693 REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: gtend,tend !< tendency … … 5644 5695 REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: root_extr_green 5645 5696 5646 REAL(wp), DIMENSION(nzb_wall:nzt_wall+1) :: lambda_green_temp !< temp. lambda5647 REAL(wp), DIMENSION(nzb_wall:nzt_wall+1) :: gamma_green_temp !< temp. gamma5697 REAL(wp), DIMENSION(nzb_wall:nzt_wall+1) :: lambda_green_temp !< temp. lambda 5698 REAL(wp), DIMENSION(nzb_wall:nzt_wall+1) :: gamma_green_temp !< temp. gamma 5648 5699 5649 5700 LOGICAL :: conserve_water_content = .true. … … 5656 5707 DO m = 1, surf_usm_h%ns 5657 5708 5658 if (surf_usm_h%frac(ind_pav_green,m).gt.0.0_wp) then 5659 ! 5660 !-- Obtain indices 5661 i = surf_usm_h%i(m) 5662 j = surf_usm_h%j(m) 5663 k = surf_usm_h%k(m) 5664 5665 DO kw = nzb_wall, nzt_wall 5666 ! 5667 !-- Calculate volumetric heat capacity of the soil, taking 5668 !-- into account water content 5669 surf_usm_h%rho_c_total_green(kw,m) = (surf_usm_h%rho_c_green(kw,m) * (1.0_wp - swc_sat_h(kw,m)) & 5670 + rho_c_water * swc_h(kw,m)) 5709 IF (surf_usm_h%frac(ind_pav_green,m) > 0.0_wp) THEN 5710 ! 5711 !-- Obtain indices 5712 i = surf_usm_h%i(m) 5713 j = surf_usm_h%j(m) 5714 k = surf_usm_h%k(m) 5671 5715 5672 ! 5673 !-- Calculate soil heat conductivity at the center of the soil 5674 !-- layers 5675 lambda_h_green_sat = lambda_h_green_sm ** (1.0_wp - swc_sat_h(kw,m)) * & 5676 lambda_h_water ** swc_h(kw,m) 5716 DO kw = nzb_wall, nzt_wall 5717 ! 5718 !-- Calculate volumetric heat capacity of the soil, taking 5719 !-- into account water content 5720 surf_usm_h%rho_c_total_green(kw,m) = (surf_usm_h%rho_c_green(kw,m) * (1.0_wp - swc_sat_h(kw,m)) & 5721 + rho_c_water * swc_h(kw,m)) 5722 5723 ! 5724 !-- Calculate soil heat conductivity at the center of the soil 5725 !-- layers 5726 lambda_h_green_sat = lambda_h_green_sm ** (1.0_wp - swc_sat_h(kw,m)) * & 5727 lambda_h_water ** swc_h(kw,m) 5728 5729 ke = 1.0_wp + LOG10(MAX(0.1_wp,swc_h(kw,m) & 5730 / swc_sat_h(kw,m))) 5731 5732 lambda_green_temp(kw) = ke * (lambda_h_green_sat - lambda_h_green_dry) + & 5733 lambda_h_green_dry 5677 5734 5678 ke = 1.0_wp + LOG10(MAX(0.1_wp,swc_h(kw,m) & 5679 / swc_sat_h(kw,m))) 5735 ENDDO 5680 5736 5681 lambda_green_temp(kw) = ke * (lambda_h_green_sat - lambda_h_green_dry) + & 5682 lambda_h_green_dry 5683 5684 ENDDO 5685 5686 5687 ! 5688 !-- Calculate soil heat conductivity (lambda_h) at the _stag level 5689 !-- using linear interpolation. For pavement surface, the 5690 !-- true pavement depth is considered 5691 DO kw = nzb_wall, nzt_wall 5692 surf_usm_h%lambda_h_green(kw,m) = ( lambda_green_temp(kw+1) + lambda_green_temp(kw) ) & 5693 * 0.5_wp 5694 ENDDO 5695 ! surf_usm_h%lambda_h_green(nzt_wall+1,m) = lambda_green_temp(nzt_wall+1) 5696 !-------------------------------------------------------------------------- 5697 5698 t_green_h(nzt_wall+1,m) = t_wall_h(nzb_wall,m) 5737 5738 ! 5739 !-- Calculate soil heat conductivity (lambda_h) at the _stag level 5740 !-- using linear interpolation. For pavement surface, the 5741 !-- true pavement depth is considered 5742 DO kw = nzb_wall, nzt_wall 5743 surf_usm_h%lambda_h_green(kw,m) = ( lambda_green_temp(kw+1) + lambda_green_temp(kw) ) & 5744 * 0.5_wp 5745 ENDDO 5746 5747 t_green_h(nzt_wall+1,m) = t_wall_h(nzb_wall,m) 5699 5748 ! 5700 5749 !-- prognostic equation for ground/roof temperature t_green_h 5701 gtend(:) = 0.0_wp5702 gtend(nzb_wall) = (1.0_wp / surf_usm_h%rho_c_total_green(nzb_wall,m)) * &5703 ( surf_usm_h%lambda_h_green(nzb_wall,m) * &5704 ( t_green_h(nzb_wall+1,m) &5705 - t_green_h(nzb_wall,m) ) * &5706 surf_usm_h%ddz_green(nzb_wall+1,m) &5707 + surf_usm_h%wghf_eb_green(m) ) * &5708 surf_usm_h%ddz_green_stag(nzb_wall,m)5709 5710 DO kw = nzb_wall+1, nzt_wall5711 gtend(kw) = (1.0_wp / surf_usm_h%rho_c_total_green(kw,m))&5712 * ( surf_usm_h%lambda_h_green(kw,m) &5713 * ( t_green_h(kw+1,m) - t_green_h(kw,m) ) &5714 * surf_usm_h%ddz_green(kw+1,m) &5715 - surf_usm_h%lambda_h_green(kw-1,m) &5716 * ( t_green_h(kw,m) - t_green_h(kw-1,m) ) &5717 * surf_usm_h%ddz_green(kw,m) &5718 ) * surf_usm_h%ddz_green_stag(kw,m)5719 ENDDO5720 5721 t_green_h_p(nzb_wall:nzt_wall,m) = t_green_h(nzb_wall:nzt_wall,m) &5722 + dt_3d * ( tsc(2) &5723 * gtend(nzb_wall:nzt_wall) + tsc(3) &5724 * surf_usm_h%tt_green_m(nzb_wall:nzt_wall,m) )5725 5726 5750 gtend(:) = 0.0_wp 5751 gtend(nzb_wall) = (1.0_wp / surf_usm_h%rho_c_total_green(nzb_wall,m)) * & 5752 ( surf_usm_h%lambda_h_green(nzb_wall,m) * & 5753 ( t_green_h(nzb_wall+1,m) & 5754 - t_green_h(nzb_wall,m) ) * & 5755 surf_usm_h%ddz_green(nzb_wall+1,m) & 5756 + surf_usm_h%wghf_eb_green(m) ) * & 5757 surf_usm_h%ddz_green_stag(nzb_wall,m) 5758 5759 DO kw = nzb_wall+1, nzt_wall 5760 gtend(kw) = (1.0_wp / surf_usm_h%rho_c_total_green(kw,m)) & 5761 * ( surf_usm_h%lambda_h_green(kw,m) & 5762 * ( t_green_h(kw+1,m) - t_green_h(kw,m) ) & 5763 * surf_usm_h%ddz_green(kw+1,m) & 5764 - surf_usm_h%lambda_h_green(kw-1,m) & 5765 * ( t_green_h(kw,m) - t_green_h(kw-1,m) ) & 5766 * surf_usm_h%ddz_green(kw,m) & 5767 ) * surf_usm_h%ddz_green_stag(kw,m) 5768 ENDDO 5769 5770 t_green_h_p(nzb_wall:nzt_wall,m) = t_green_h(nzb_wall:nzt_wall,m) & 5771 + dt_3d * ( tsc(2) & 5772 * gtend(nzb_wall:nzt_wall) + tsc(3) & 5773 * surf_usm_h%tt_green_m(nzb_wall:nzt_wall,m) ) 5774 5775 5727 5776 ! 5728 5777 !-- calculate t_green tendencies for the next Runge-Kutta step 5729 IF ( timestep_scheme(1:5) == 'runge' ) THEN5730 IF ( intermediate_timestep_count == 1 ) THEN5731 DO kw = nzb_wall, nzt_wall5732 surf_usm_h%tt_green_m(kw,m) = gtend(kw)5733 ENDDO5734 ELSEIF ( intermediate_timestep_count < &5735 intermediate_timestep_count_max ) THEN5736 DO kw = nzb_wall, nzt_wall5737 surf_usm_h%tt_green_m(kw,m) = -9.5625_wp * gtend(kw) + &5738 5.3125_wp * surf_usm_h%tt_green_m(kw,m)5739 ENDDO5740 ENDIF5741 ENDIF5742 5743 !--------------------------------------------------------------5744 DO kw = nzb_wall, nzt_wall5745 5746 !5747 !-- Calculate soil diffusivity at the center of the soil layers5748 lambda_green_temp(kw) = (- b_ch * surf_usm_h%gamma_w_green_sat(kw,m) * psi_sat &5749 / swc_sat_h(kw,m) ) * ( MAX( swc_h(kw,m), &5750 wilt_h(kw,m) ) / swc_sat_h(kw,m) )**( &5751 b_ch + 2.0_wp )5752 5753 !5754 !-- Parametrization of Van Genuchten5755 IF ( soil_type /= 7 ) THEN5756 !5757 !-- Calculate the hydraulic conductivity after Van Genuchten5758 !-- (1980)5759 h_vg = ( ( (swc_res_h(kw,m) - swc_sat_h(kw,m)) / ( swc_res_h(kw,m) - &5760 MAX( swc_h(kw,m), wilt_h(kw,m) ) ) )**( &5761 surf_usm_h%n_vg_green(m) / (surf_usm_h%n_vg_green(m) - 1.0_wp ) ) - 1.0_wp &5762 )**( 1.0_wp / surf_usm_h%n_vg_green(m) ) / surf_usm_h%alpha_vg_green(m)5763 5764 5765 gamma_green_temp(kw) = surf_usm_h%gamma_w_green_sat(kw,m) * ( ( (1.0_wp + &5766 ( surf_usm_h%alpha_vg_green(m) * h_vg )**surf_usm_h%n_vg_green(m))**( &5767 1.0_wp - 1.0_wp / surf_usm_h%n_vg_green(m) ) - ( &5768 surf_usm_h%alpha_vg_green(m) * h_vg )**( surf_usm_h%n_vg_green(m) &5769 - 1.0_wp) )**2 ) &5770 / ( ( 1.0_wp + ( surf_usm_h%alpha_vg_green(m) * h_vg &5771 )**surf_usm_h%n_vg_green(m) )**( ( 1.0_wp - 1.0_wp &5772 / surf_usm_h%n_vg_green(m) ) *( surf_usm_h%l_vg_green(m) + 2.0_wp) ) )5773 5774 !5775 !-- Parametrization of Clapp & Hornberger5776 ELSE5777 gamma_green_temp(kw) = surf_usm_h%gamma_w_green_sat(kw,m) * ( swc_h(kw,m) &5778 / swc_sat_h(kw,m) )**(2.0_wp * b_ch + 3.0_wp)5779 ENDIF5780 5781 ENDDO5782 5783 !5784 !-- Prognostic equation for soil moisture content. Only performed,5785 !-- when humidity is enabled in the atmosphere5786 IF ( humidity ) THEN5787 !5788 !-- Calculate soil diffusivity (lambda_w) at the _stag level5789 !-- using linear interpolation. To do: replace this with5790 !-- ECMWF-IFS Eq. 8.815791 DO kw = nzb_wall, nzt_wall-15792 5793 surf_usm_h%lambda_w_green(kw,m) = ( lambda_green_temp(kw+1) + lambda_green_temp(kw) ) &5794 * 0.5_wp5795 surf_usm_h%gamma_w_green(kw,m) = ( gamma_green_temp(kw+1) + gamma_green_temp(kw) ) &5796 * 0.5_wp5797 5798 ENDDO5799 5800 !5801 !5802 !-- In case of a closed bottom (= water content is conserved),5803 !-- set hydraulic conductivity to zero to that no water will be5804 !-- lost in the bottom layer.5805 IF ( conserve_water_content ) THEN5806 surf_usm_h%gamma_w_green(kw,m) = 0.0_wp5807 ELSE5808 surf_usm_h%gamma_w_green(kw,m) = gamma_green_temp(nzt_wall)5809 ENDIF5810 5811 !-- The root extraction (= root_extr * qsws_veg / (rho_l5812 !-- * l_v)) ensures the mass conservation for water. The5813 !-- transpiration of plants equals the cumulative withdrawals by5814 !-- the roots in the soil. The scheme takes into account the5815 !-- availability of water in the soil layers as well as the root5816 !-- fraction in the respective layer. Layer with moisture below5817 !-- wilting point will not contribute, which reflects the5818 !-- preference of plants to take water from moister layers.5819 5820 !5821 !-- Calculate the root extraction (ECMWF 7.69, the sum of5822 !-- root_extr = 1). The energy balance solver guarantees a5823 !-- positive transpiration, so that there is no need for an5824 !-- additional check.5825 m_total = 0.0_wp5826 DO kw = nzb_wall, nzt_wall5827 IF ( swc_h(kw,m) > wilt_h(kw,m) ) THEN5828 m_total = m_total + rootfr_h(kw,m) * swc_h(kw,m)5829 ENDIF5830 ENDDO5831 5832 IF ( m_total > 0.0_wp ) THEN5833 DO kw = nzb_wall, nzt_wall5834 IF ( swc_h(kw,m) > wilt_h(kw,m) ) THEN5835 root_extr_green(kw) = rootfr_h(kw,m) * swc_h(kw,m) &5836 / m_total5837 ELSE5838 root_extr_green(kw) = 0.0_wp5839 ENDIF5840 ENDDO5841 ENDIF5842 5843 !5844 !-- Prognostic equation for soil water content m_soil.5845 tend(:) = 0.0_wp5846 5847 tend(nzb_wall) = ( surf_usm_h%lambda_w_green(nzb_wall,m) * ( &5848 swc_h(nzb_wall+1,m) - swc_h(nzb_wall,m) ) &5849 * surf_usm_h%ddz_green(nzb_wall+1,m) - surf_usm_h%gamma_w_green(nzb_wall,m) - ( &5850 root_extr_green(nzb_wall) * surf_usm_h%qsws_veg(m) &5851 ! + surf_usm_h%qsws_soil_green(m)5852 ) * drho_l_lv ) &5853 * surf_usm_h%ddz_green_stag(nzb_wall,m)5854 5855 DO kw = nzb_wall+1, nzt_wall-15856 tend(kw) = ( surf_usm_h%lambda_w_green(kw,m) * ( swc_h(kw+1,m) &5857 - swc_h(kw,m) ) * surf_usm_h%ddz_green(kw+1,m) &5858 - surf_usm_h%gamma_w_green(kw,m) &5859 - surf_usm_h%lambda_w_green(kw-1,m) * (swc_h(kw,m) - &5860 swc_h(kw-1,m)) * surf_usm_h%ddz_green(kw,m) &5861 + surf_usm_h%gamma_w_green(kw-1,m) - (root_extr_green(kw) &5862 * surf_usm_h%qsws_veg(m) * drho_l_lv) &5863 ) * surf_usm_h%ddz_green_stag(kw,m)5864 5865 ENDDO5866 tend(nzt_wall) = ( - surf_usm_h%gamma_w_green(nzt_wall,m) &5867 - surf_usm_h%lambda_w_green(nzt_wall-1,m) &5868 * (swc_h(nzt_wall,m) &5869 - swc_h(nzt_wall-1,m)) &5870 * surf_usm_h%ddz_green(nzt_wall,m) &5871 + surf_usm_h%gamma_w_green(nzt_wall-1,m) - ( &5872 root_extr_green(nzt_wall) &5873 * surf_usm_h%qsws_veg(m) * drho_l_lv ) &5874 ) * surf_usm_h%ddz_green_stag(nzt_wall,m)5875 5876 swc_h_p(nzb_wall:nzt_wall,m) = swc_h(nzb_wall:nzt_wall,m)&5877 + dt_3d * ( tsc(2) * tend(:) &5878 + tsc(3) * surf_usm_h%tswc_h_m(:,m) )5879 5880 !5881 !-- Account for dry soils (find a better solution here!)5882 DO kw = nzb_wall, nzt_wall5883 IF ( swc_h_p(kw,m) < 0.0_wp ) swc_h_p(kw,m) = 0.0_wp5884 ENDDO5885 5886 !5887 !-- Calculate m_soil tendencies for the next Runge-Kutta step5888 IF ( timestep_scheme(1:5) == 'runge' ) THEN5889 IF ( intermediate_timestep_count == 1 ) THEN5890 DO kw = nzb_wall, nzt_wall5891 surf_usm_h%tswc_h_m(kw,m) = tend(kw)5892 ENDDO5893 ELSEIF ( intermediate_timestep_count < &5894 intermediate_timestep_count_max ) THEN5895 DO kw = nzb_wall, nzt_wall5896 surf_usm_h%tswc_h_m(kw,m) = -9.5625_wp * tend(kw) + 5.3125_wp&5897 * surf_usm_h%tswc_h_m(kw,m)5898 ENDDO5899 ENDIF5900 ENDIF5901 ENDIF5902 !--------------------------------------------------------------5903 ENDIF5904 5905 ENDDO5906 5907 !5908 !-- For vertical surfaces5909 DO l = 0, 35910 DO m = 1, surf_usm_v(l)%ns5911 5912 if (surf_usm_v(l)%frac(ind_pav_green,m).gt.0.0_wp) then5913 if (1.gt.2) then5914 !5915 !-- Obtain indices5916 i = surf_usm_v(l)%i(m)5917 j = surf_usm_v(l)%j(m)5918 k = surf_usm_v(l)%k(m)5919 5920 t_green_v(l)%t(nzt_wall+1,m) = t_wall_v(l)%t(nzb_wall,m)5921 !5922 !-- prognostic equation for green temperature t_green_v5923 gtend(:) = 0.0_wp5924 gtend(nzb_wall) = (1.0_wp / surf_usm_v(l)%rho_c_green(nzb_wall,m)) * &5925 ( surf_usm_v(l)%lambda_h_green(nzb_wall,m) * &5926 ( t_green_v(l)%t(nzb_wall+1,m) &5927 - t_green_v(l)%t(nzb_wall,m) ) * &5928 surf_usm_v(l)%ddz_green(nzb_wall+1,m) &5929 + surf_usm_v(l)%wghf_eb(m) ) * &5930 surf_usm_v(l)%ddz_green_stag(nzb_wall,m)5931 5932 DO kw = nzb_wall+1, nzt_wall5933 gtend(kw) = (1.0_wp / surf_usm_v(l)%rho_c_green(kw,m)) &5934 * ( surf_usm_v(l)%lambda_h_green(kw,m) &5935 * ( t_green_v(l)%t(kw+1,m) - t_green_v(l)%t(kw,m) ) &5936 * surf_usm_v(l)%ddz_green(kw+1,m) &5937 - surf_usm_v(l)%lambda_h(kw-1,m) &5938 * ( t_green_v(l)%t(kw,m) - t_green_v(l)%t(kw-1,m) ) &5939 * surf_usm_v(l)%ddz_green(kw,m) ) &5940 * surf_usm_v(l)%ddz_green_stag(kw,m)5941 ENDDO5942 5943 t_green_v_p(l)%t(nzb_wall:nzt_wall,m) = &5944 t_green_v(l)%t(nzb_wall:nzt_wall,m) &5945 + dt_3d * ( tsc(2) &5946 * gtend(nzb_wall:nzt_wall) + tsc(3) &5947 * surf_usm_v(l)%tt_green_m(nzb_wall:nzt_wall,m) )5948 5949 !5950 !-- calculate t_green tendencies for the next Runge-Kutta step5951 5778 IF ( timestep_scheme(1:5) == 'runge' ) THEN 5952 5779 IF ( intermediate_timestep_count == 1 ) THEN 5953 5780 DO kw = nzb_wall, nzt_wall 5954 surf_usm_ v(l)%tt_green_m(kw,m) = gtend(kw)5781 surf_usm_h%tt_green_m(kw,m) = gtend(kw) 5955 5782 ENDDO 5956 5783 ELSEIF ( intermediate_timestep_count < & 5957 5784 intermediate_timestep_count_max ) THEN 5958 5785 DO kw = nzb_wall, nzt_wall 5959 surf_usm_v(l)%tt_green_m(kw,m) = & 5960 - 9.5625_wp * gtend(kw) + & 5961 5.3125_wp * surf_usm_v(l)%tt_green_m(kw,m) 5786 surf_usm_h%tt_green_m(kw,m) = -9.5625_wp * gtend(kw) + & 5787 5.3125_wp * surf_usm_h%tt_green_m(kw,m) 5962 5788 ENDDO 5963 5789 ENDIF 5964 5790 ENDIF 5965 endif 5966 5967 DO kw = nzb_wall, nzt_wall+1 5968 t_green_v(l)%t(kw,m) = t_wall_v(l)%t(nzb_wall,m) 5969 ENDDO 5791 5792 DO kw = nzb_wall, nzt_wall 5793 5794 ! 5795 !-- Calculate soil diffusivity at the center of the soil layers 5796 lambda_green_temp(kw) = (- b_ch * surf_usm_h%gamma_w_green_sat(kw,m) * psi_sat & 5797 / swc_sat_h(kw,m) ) * ( MAX( swc_h(kw,m), & 5798 wilt_h(kw,m) ) / swc_sat_h(kw,m) )**( & 5799 b_ch + 2.0_wp ) 5800 5801 ! 5802 !-- Parametrization of Van Genuchten 5803 IF ( soil_type /= 7 ) THEN 5804 ! 5805 !-- Calculate the hydraulic conductivity after Van Genuchten 5806 !-- (1980) 5807 h_vg = ( ( (swc_res_h(kw,m) - swc_sat_h(kw,m)) / ( swc_res_h(kw,m) - & 5808 MAX( swc_h(kw,m), wilt_h(kw,m) ) ) )**( & 5809 surf_usm_h%n_vg_green(m) / (surf_usm_h%n_vg_green(m) - 1.0_wp ) ) - 1.0_wp & 5810 )**( 1.0_wp / surf_usm_h%n_vg_green(m) ) / surf_usm_h%alpha_vg_green(m) 5811 5812 5813 gamma_green_temp(kw) = surf_usm_h%gamma_w_green_sat(kw,m) * ( ( (1.0_wp + & 5814 ( surf_usm_h%alpha_vg_green(m) * h_vg )**surf_usm_h%n_vg_green(m))**( & 5815 1.0_wp - 1.0_wp / surf_usm_h%n_vg_green(m) ) - ( & 5816 surf_usm_h%alpha_vg_green(m) * h_vg )**( surf_usm_h%n_vg_green(m) & 5817 - 1.0_wp) )**2 ) & 5818 / ( ( 1.0_wp + ( surf_usm_h%alpha_vg_green(m) * h_vg & 5819 )**surf_usm_h%n_vg_green(m) )**( ( 1.0_wp - 1.0_wp & 5820 / surf_usm_h%n_vg_green(m) ) *( surf_usm_h%l_vg_green(m) + 2.0_wp) ) ) 5821 5822 ! 5823 !-- Parametrization of Clapp & Hornberger 5824 ELSE 5825 gamma_green_temp(kw) = surf_usm_h%gamma_w_green_sat(kw,m) * ( swc_h(kw,m) & 5826 / swc_sat_h(kw,m) )**(2.0_wp * b_ch + 3.0_wp) 5827 ENDIF 5828 5829 ENDDO 5830 5831 ! 5832 !-- Prognostic equation for soil moisture content. Only performed, 5833 !-- when humidity is enabled in the atmosphere 5834 IF ( humidity ) THEN 5835 ! 5836 !-- Calculate soil diffusivity (lambda_w) at the _stag level 5837 !-- using linear interpolation. To do: replace this with 5838 !-- ECMWF-IFS Eq. 8.81 5839 DO kw = nzb_wall, nzt_wall-1 5840 5841 surf_usm_h%lambda_w_green(kw,m) = ( lambda_green_temp(kw+1) + lambda_green_temp(kw) ) & 5842 * 0.5_wp 5843 surf_usm_h%gamma_w_green(kw,m) = ( gamma_green_temp(kw+1) + gamma_green_temp(kw) ) & 5844 * 0.5_wp 5845 5846 ENDDO 5847 5848 ! 5849 !-- In case of a closed bottom (= water content is conserved), 5850 !-- set hydraulic conductivity to zero to that no water will be 5851 !-- lost in the bottom layer. 5852 IF ( conserve_water_content ) THEN 5853 surf_usm_h%gamma_w_green(kw,m) = 0.0_wp 5854 ELSE 5855 surf_usm_h%gamma_w_green(kw,m) = gamma_green_temp(nzt_wall) 5856 ENDIF 5857 5858 !-- The root extraction (= root_extr * qsws_veg / (rho_l 5859 !-- * l_v)) ensures the mass conservation for water. The 5860 !-- transpiration of plants equals the cumulative withdrawals by 5861 !-- the roots in the soil. The scheme takes into account the 5862 !-- availability of water in the soil layers as well as the root 5863 !-- fraction in the respective layer. Layer with moisture below 5864 !-- wilting point will not contribute, which reflects the 5865 !-- preference of plants to take water from moister layers. 5866 5867 ! 5868 !-- Calculate the root extraction (ECMWF 7.69, the sum of 5869 !-- root_extr = 1). The energy balance solver guarantees a 5870 !-- positive transpiration, so that there is no need for an 5871 !-- additional check. 5872 m_total = 0.0_wp 5873 DO kw = nzb_wall, nzt_wall 5874 IF ( swc_h(kw,m) > wilt_h(kw,m) ) THEN 5875 m_total = m_total + rootfr_h(kw,m) * swc_h(kw,m) 5876 ENDIF 5877 ENDDO 5878 5879 IF ( m_total > 0.0_wp ) THEN 5880 DO kw = nzb_wall, nzt_wall 5881 IF ( swc_h(kw,m) > wilt_h(kw,m) ) THEN 5882 root_extr_green(kw) = rootfr_h(kw,m) * swc_h(kw,m) & 5883 / m_total 5884 ELSE 5885 root_extr_green(kw) = 0.0_wp 5886 ENDIF 5887 ENDDO 5888 ENDIF 5889 5890 ! 5891 !-- Prognostic equation for soil water content m_soil. 5892 tend(:) = 0.0_wp 5893 5894 tend(nzb_wall) = ( surf_usm_h%lambda_w_green(nzb_wall,m) * ( & 5895 swc_h(nzb_wall+1,m) - swc_h(nzb_wall,m) ) & 5896 * surf_usm_h%ddz_green(nzb_wall+1,m) - surf_usm_h%gamma_w_green(nzb_wall,m) - ( & 5897 root_extr_green(nzb_wall) * surf_usm_h%qsws_veg(m) & 5898 ! + surf_usm_h%qsws_soil_green(m) 5899 ) * drho_l_lv ) & 5900 * surf_usm_h%ddz_green_stag(nzb_wall,m) 5901 5902 DO kw = nzb_wall+1, nzt_wall-1 5903 tend(kw) = ( surf_usm_h%lambda_w_green(kw,m) * ( swc_h(kw+1,m) & 5904 - swc_h(kw,m) ) * surf_usm_h%ddz_green(kw+1,m) & 5905 - surf_usm_h%gamma_w_green(kw,m) & 5906 - surf_usm_h%lambda_w_green(kw-1,m) * (swc_h(kw,m) - & 5907 swc_h(kw-1,m)) * surf_usm_h%ddz_green(kw,m) & 5908 + surf_usm_h%gamma_w_green(kw-1,m) - (root_extr_green(kw) & 5909 * surf_usm_h%qsws_veg(m) * drho_l_lv) & 5910 ) * surf_usm_h%ddz_green_stag(kw,m) 5911 5912 ENDDO 5913 tend(nzt_wall) = ( - surf_usm_h%gamma_w_green(nzt_wall,m) & 5914 - surf_usm_h%lambda_w_green(nzt_wall-1,m) & 5915 * (swc_h(nzt_wall,m) & 5916 - swc_h(nzt_wall-1,m)) & 5917 * surf_usm_h%ddz_green(nzt_wall,m) & 5918 + surf_usm_h%gamma_w_green(nzt_wall-1,m) - ( & 5919 root_extr_green(nzt_wall) & 5920 * surf_usm_h%qsws_veg(m) * drho_l_lv ) & 5921 ) * surf_usm_h%ddz_green_stag(nzt_wall,m) 5922 5923 swc_h_p(nzb_wall:nzt_wall,m) = swc_h(nzb_wall:nzt_wall,m)& 5924 + dt_3d * ( tsc(2) * tend(:) & 5925 + tsc(3) * surf_usm_h%tswc_h_m(:,m) ) 5926 5927 ! 5928 !-- Account for dry soils (find a better solution here!) 5929 DO kw = nzb_wall, nzt_wall 5930 IF ( swc_h_p(kw,m) < 0.0_wp ) swc_h_p(kw,m) = 0.0_wp 5931 ENDDO 5932 5933 ! 5934 !-- Calculate m_soil tendencies for the next Runge-Kutta step 5935 IF ( timestep_scheme(1:5) == 'runge' ) THEN 5936 IF ( intermediate_timestep_count == 1 ) THEN 5937 DO kw = nzb_wall, nzt_wall 5938 surf_usm_h%tswc_h_m(kw,m) = tend(kw) 5939 ENDDO 5940 ELSEIF ( intermediate_timestep_count < & 5941 intermediate_timestep_count_max ) THEN 5942 DO kw = nzb_wall, nzt_wall 5943 surf_usm_h%tswc_h_m(kw,m) = -9.5625_wp * tend(kw) + 5.3125_wp& 5944 * surf_usm_h%tswc_h_m(kw,m) 5945 ENDDO 5946 ENDIF 5947 ENDIF 5948 ENDIF 5949 5950 ENDIF 5951 5952 ENDDO 5953 5954 ! 5955 !-- For vertical surfaces 5956 DO l = 0, 3 5957 DO m = 1, surf_usm_v(l)%ns 5958 5959 IF (surf_usm_v(l)%frac(ind_pav_green,m) > 0.0_wp) THEN 5960 ! 5961 !-- no substrate layer for green walls / only groundbase green walls (ivy i.e.) -> green layers get same 5962 !-- temperature as first wall layer 5963 !-- there fore no temperature calculations for vertical green substrate layers now 5964 5965 ! 5966 ! ! 5967 ! !-- Obtain indices 5968 ! i = surf_usm_v(l)%i(m) 5969 ! j = surf_usm_v(l)%j(m) 5970 ! k = surf_usm_v(l)%k(m) 5971 ! 5972 ! t_green_v(l)%t(nzt_wall+1,m) = t_wall_v(l)%t(nzb_wall,m) 5973 ! ! 5974 ! !-- prognostic equation for green temperature t_green_v 5975 ! gtend(:) = 0.0_wp 5976 ! gtend(nzb_wall) = (1.0_wp / surf_usm_v(l)%rho_c_green(nzb_wall,m)) * & 5977 ! ( surf_usm_v(l)%lambda_h_green(nzb_wall,m) * & 5978 ! ( t_green_v(l)%t(nzb_wall+1,m) & 5979 ! - t_green_v(l)%t(nzb_wall,m) ) * & 5980 ! surf_usm_v(l)%ddz_green(nzb_wall+1,m) & 5981 ! + surf_usm_v(l)%wghf_eb(m) ) * & 5982 ! surf_usm_v(l)%ddz_green_stag(nzb_wall,m) 5983 ! 5984 ! DO kw = nzb_wall+1, nzt_wall 5985 ! gtend(kw) = (1.0_wp / surf_usm_v(l)%rho_c_green(kw,m)) & 5986 ! * ( surf_usm_v(l)%lambda_h_green(kw,m) & 5987 ! * ( t_green_v(l)%t(kw+1,m) - t_green_v(l)%t(kw,m) ) & 5988 ! * surf_usm_v(l)%ddz_green(kw+1,m) & 5989 ! - surf_usm_v(l)%lambda_h(kw-1,m) & 5990 ! * ( t_green_v(l)%t(kw,m) - t_green_v(l)%t(kw-1,m) ) & 5991 ! * surf_usm_v(l)%ddz_green(kw,m) ) & 5992 ! * surf_usm_v(l)%ddz_green_stag(kw,m) 5993 ! ENDDO 5994 ! 5995 ! t_green_v_p(l)%t(nzb_wall:nzt_wall,m) = & 5996 ! t_green_v(l)%t(nzb_wall:nzt_wall,m) & 5997 ! + dt_3d * ( tsc(2) & 5998 ! * gtend(nzb_wall:nzt_wall) + tsc(3) & 5999 ! * surf_usm_v(l)%tt_green_m(nzb_wall:nzt_wall,m) ) 6000 ! 6001 ! ! 6002 ! !-- calculate t_green tendencies for the next Runge-Kutta step 6003 ! IF ( timestep_scheme(1:5) == 'runge' ) THEN 6004 ! IF ( intermediate_timestep_count == 1 ) THEN 6005 ! DO kw = nzb_wall, nzt_wall 6006 ! surf_usm_v(l)%tt_green_m(kw,m) = gtend(kw) 6007 ! ENDDO 6008 ! ELSEIF ( intermediate_timestep_count < & 6009 ! intermediate_timestep_count_max ) THEN 6010 ! DO kw = nzb_wall, nzt_wall 6011 ! surf_usm_v(l)%tt_green_m(kw,m) = & 6012 ! - 9.5625_wp * gtend(kw) + & 6013 ! 5.3125_wp * surf_usm_v(l)%tt_green_m(kw,m) 6014 ! ENDDO 6015 ! ENDIF 6016 ! ENDIF 6017 6018 DO kw = nzb_wall, nzt_wall+1 6019 t_green_v(l)%t(kw,m) = t_wall_v(l)%t(nzb_wall,m) 6020 ENDDO 5970 6021 5971 ENDIF6022 ENDIF 5972 6023 5973 6024 ENDDO … … 6025 6076 usm_wall_mod 6026 6077 6027 6028 6078 6029 6079 ! … … 6094 6144 SUBROUTINE usm_read_anthropogenic_heat 6095 6145 6096 INTEGER(iwp) :: i,j,k,ii 6097 REAL(wp) :: heat 6098 6146 INTEGER(iwp) :: i,j,k,ii !< running indices 6147 REAL(wp) :: heat !< anthropogenic heat 6148 6149 ! 6099 6150 !-- allocation of array of sources of anthropogenic heat and their diural profile 6100 6151 ALLOCATE( aheat(naheatlayers,nys:nyn,nxl:nxr) ) 6101 6152 ALLOCATE( aheatprof(naheatlayers,0:24) ) 6102 6153 6103 ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!6154 ! 6104 6155 !-- read daily amount of heat and its daily cycle 6105 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!6106 6156 aheat = 0.0_wp 6107 6157 DO ii = 0, io_blocks-1 … … 6136 6186 ENDDO 6137 6187 6138 ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!6188 ! 6139 6189 !-- read diurnal profiles of heat sources 6140 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!6141 6190 aheatprof = 0.0_wp 6142 6191 DO ii = 0, io_blocks-1 6143 6192 IF ( ii == io_group ) THEN 6144 6193 ! 6145 6194 !-- open anthropogenic heat profile file 6146 6195 OPEN( 151, file='ANTHROPOGENIC_HEAT_PROFILE'//TRIM(coupling_char), action='read', & … … 6188 6237 IMPLICIT NONE 6189 6238 6190 INTEGER(iwp) :: l !< index variable for surface type6191 INTEGER(iwp) :: i !< running index over input files6192 INTEGER(iwp) :: k !< running index over previous input files covering current local domain6193 INTEGER(iwp) :: ns_h_on_file_usm !< number of horizontal surface elements (urban type) on file6194 INTEGER(iwp) :: nxlc !< index of left boundary on current subdomain6195 INTEGER(iwp) :: nxlf !< index of left boundary on former subdomain6196 INTEGER(iwp) :: nxl_on_file !< index of left boundary on former local domain6197 INTEGER(iwp) :: nxrc !< index of right boundary on current subdomain6198 INTEGER(iwp) :: nxrf !< index of right boundary on former subdomain6199 INTEGER(iwp) :: nxr_on_file !< index of right boundary on former local domain6200 INTEGER(iwp) :: nync !< index of north boundary on current subdomain6201 INTEGER(iwp) :: nynf !< index of north boundary on former subdomain6202 INTEGER(iwp) :: nyn_on_file !< index of north boundary on former local domain6203 INTEGER(iwp) :: nysc !< index of south boundary on current subdomain6204 INTEGER(iwp) :: nysf !< index of south boundary on former subdomain6205 INTEGER(iwp) :: nys_on_file !< index of south boundary on former local domain6239 INTEGER(iwp) :: l !< index variable for surface type 6240 INTEGER(iwp) :: i !< running index over input files 6241 INTEGER(iwp) :: k !< running index over previous input files covering current local domain 6242 INTEGER(iwp) :: ns_h_on_file_usm !< number of horizontal surface elements (urban type) on file 6243 INTEGER(iwp) :: nxlc !< index of left boundary on current subdomain 6244 INTEGER(iwp) :: nxlf !< index of left boundary on former subdomain 6245 INTEGER(iwp) :: nxl_on_file !< index of left boundary on former local domain 6246 INTEGER(iwp) :: nxrc !< index of right boundary on current subdomain 6247 INTEGER(iwp) :: nxrf !< index of right boundary on former subdomain 6248 INTEGER(iwp) :: nxr_on_file !< index of right boundary on former local domain 6249 INTEGER(iwp) :: nync !< index of north boundary on current subdomain 6250 INTEGER(iwp) :: nynf !< index of north boundary on former subdomain 6251 INTEGER(iwp) :: nyn_on_file !< index of north boundary on former local domain 6252 INTEGER(iwp) :: nysc !< index of south boundary on current subdomain 6253 INTEGER(iwp) :: nysf !< index of south boundary on former subdomain 6254 INTEGER(iwp) :: nys_on_file !< index of south boundary on former local domain 6206 6255 6207 INTEGER(iwp) :: ns_v_on_file_usm(0:3) !< number of vertical surface elements (urban type) on file6256 INTEGER(iwp) :: ns_v_on_file_usm(0:3) !< number of vertical surface elements (urban type) on file 6208 6257 6209 6258 INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE, SAVE :: start_index_on_file … … 6260 6309 6261 6310 DO l = 0, 3 6262 IF ( ALLOCATED( tmp_surf_wall_v(l)%t ) ) 6311 IF ( ALLOCATED( tmp_surf_wall_v(l)%t ) ) & 6263 6312 DEALLOCATE( tmp_surf_wall_v(l)%t ) 6264 6313 IF ( ALLOCATED( tmp_wall_v(l)%t ) ) & … … 6322 6371 CASE ( 't_surf_wall_h' ) 6323 6372 IF ( k == 1 ) THEN 6324 IF ( .NOT. ALLOCATED( t_surf_wall_h_1 ) ) 6373 IF ( .NOT. ALLOCATED( t_surf_wall_h_1 ) ) & 6325 6374 ALLOCATE( t_surf_wall_h_1(1:surf_usm_h%ns) ) 6326 6375 READ ( 13 ) tmp_surf_wall_h … … 6338 6387 CASE ( 't_surf_wall_v(0)' ) 6339 6388 IF ( k == 1 ) THEN 6340 IF ( .NOT. ALLOCATED( t_surf_wall_v_1(0)%t ) ) 6389 IF ( .NOT. ALLOCATED( t_surf_wall_v_1(0)%t ) ) & 6341 6390 ALLOCATE( t_surf_wall_v_1(0)%t(1:surf_usm_v(0)%ns) ) 6342 6391 READ ( 13 ) tmp_surf_wall_v(0)%t … … 6354 6403 CASE ( 't_surf_wall_v(1)' ) 6355 6404 IF ( k == 1 ) THEN 6356 IF ( .NOT. ALLOCATED( t_surf_wall_v_1(1)%t ) ) 6405 IF ( .NOT. ALLOCATED( t_surf_wall_v_1(1)%t ) ) & 6357 6406 ALLOCATE( t_surf_wall_v_1(1)%t(1:surf_usm_v(1)%ns) ) 6358 6407 READ ( 13 ) tmp_surf_wall_v(1)%t … … 6370 6419 CASE ( 't_surf_wall_v(2)' ) 6371 6420 IF ( k == 1 ) THEN 6372 IF ( .NOT. ALLOCATED( t_surf_wall_v_1(2)%t ) ) 6421 IF ( .NOT. ALLOCATED( t_surf_wall_v_1(2)%t ) ) & 6373 6422 ALLOCATE( t_surf_wall_v_1(2)%t(1:surf_usm_v(2)%ns) ) 6374 6423 READ ( 13 ) tmp_surf_wall_v(2)%t … … 6386 6435 CASE ( 't_surf_wall_v(3)' ) 6387 6436 IF ( k == 1 ) THEN 6388 IF ( .NOT. ALLOCATED( t_surf_wall_v_1(3)%t ) ) 6437 IF ( .NOT. ALLOCATED( t_surf_wall_v_1(3)%t ) ) & 6389 6438 ALLOCATE( t_surf_wall_v_1(3)%t(1:surf_usm_v(3)%ns) ) 6390 6439 READ ( 13 ) tmp_surf_wall_v(3)%t … … 6835 6884 6836 6885 END SUBROUTINE usm_rrd_local 6837 6838 6886 6839 6887 … … 6883 6931 IF ( .NOT. ascii_file ) RETURN 6884 6932 6885 ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!6933 ! 6886 6934 !-- read categories of walls and their parameters 6887 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!6888 6935 DO ii = 0, io_blocks-1 6889 6936 IF ( ii == io_group ) THEN 6890 6937 ! 6891 6938 !-- open urban surface file 6892 6939 OPEN( 151, file='SURFACE_PARAMETERS'//coupling_char, action='read', & 6893 status='old', form='formatted', err=15 ) 6940 status='old', form='formatted', err=15 ) 6941 ! 6894 6942 !-- first test and get n_surface_types 6895 6943 k = 0 … … 6906 6954 ALLOCATE( surface_type_codes(n_surface_types) ) 6907 6955 ALLOCATE( surface_params(n_surface_params, n_surface_types) ) 6956 ! 6908 6957 !-- real reading 6909 6958 rewind( 151 ) … … 6927 6976 ENDDO 6928 6977 6929 ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!6978 ! 6930 6979 !-- read types of surfaces 6931 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!6932 6980 usm_par = 0 6933 6981 DO ii = 0, io_blocks-1 6934 6982 IF ( ii == io_group ) THEN 6935 6983 6936 6984 ! 6937 6985 !-- open csv urban surface file 6938 6986 OPEN( 151, file='URBAN_SURFACE'//TRIM(coupling_char), action='read', & … … 6942 6990 DO 6943 6991 l = l+1 6992 ! 6944 6993 !-- i, j, height, nz, roof, dirwe, dirsn, category, soilcat, 6945 6994 !-- weheight1, wecat1, snheight1, sncat1, weheight2, wecat2, snheight2, sncat2, … … 6955 7004 6956 7005 IF ( i >= nxlg .AND. i <= nxrg .AND. j >= nysg .AND. j <= nyng ) THEN 7006 ! 6957 7007 !-- write integer variables into array 6958 7008 usm_par(:,j,i) = (/1, nz, roof, dirwe, dirsn, category, & 6959 7009 weheight1, wecat1, weheight2, wecat2, weheight3, wecat3, & 6960 7010 snheight1, sncat1, snheight2, sncat2, snheight3, sncat3 /) 7011 ! 6961 7012 !-- write real values into array 6962 7013 usm_val(:,j,i) = (/ albedo, thick, & … … 7003 7054 usm_par(16,j,i) > nzt & !< incorrect wall or roof level height for south-north wall 7004 7055 ) ) THEN 7056 ! 7005 7057 !-- incorrect input data 7006 7058 WRITE (message_string, "(A,2I5)") 'missing or incorrect data in file URBAN_SURFACE'// & … … 7063 7115 ENDDO 7064 7116 IF ( ip == -99999 ) THEN 7117 ! 7065 7118 !-- land/roof category not found 7066 7119 WRITE (9,"(A,I5,A,3I5)") 'land/roof category ', it, & … … 7079 7132 ENDDO 7080 7133 IF ( ip == -99999 ) THEN 7134 ! 7081 7135 !-- default land/roof category not found 7082 7136 WRITE (9,"(A,I5,A,3I5)") 'Default land/roof category', category, ' not found!' … … 7090 7144 surf_usm_h%albedo(:,m) = surface_params(ialbedo,ip) 7091 7145 ENDIF 7146 ! 7092 7147 !-- Albedo type is 0 (custom), others are replaced later 7093 7148 surf_usm_h%albedo_type(:,m) = 0 7149 ! 7094 7150 !-- Transmissivity 7095 7151 IF ( surf_usm_h%transmissivity(m) < 0.0_wp ) THEN … … 7175 7231 7176 7232 IF ( iw < 0 .OR. jw < 0 ) THEN 7233 ! 7177 7234 !-- wall on west or south border of the domain - assign default category 7178 7235 IF ( kw <= roof_height_limit ) THEN … … 7181 7238 surf_usm_v(l)%surface_types(m) = roof_category !< default category for wall surface in roof zone 7182 7239 END IF 7183 surf_usm_v(l)%albedo(:,m) = -1.0_wp7184 surf_usm_v(l)%thickness_wall(m) = -1.0_wp7185 surf_usm_v(l)%thickness_window(m) 7186 surf_usm_v(l)%thickness_green(m) 7187 surf_usm_v(l)%transmissivity(m) = -1.0_wp7240 surf_usm_v(l)%albedo(:,m) = -1.0_wp 7241 surf_usm_v(l)%thickness_wall(m) = -1.0_wp 7242 surf_usm_v(l)%thickness_window(m) = -1.0_wp 7243 surf_usm_v(l)%thickness_green(m) = -1.0_wp 7244 surf_usm_v(l)%transmissivity(m) = -1.0_wp 7188 7245 ELSE IF ( kw <= usm_par(ii,jw,iw) ) THEN 7246 ! 7189 7247 !-- pedestrian zone 7190 7248 IF ( usm_par(ii+1,jw,iw) == 0 ) THEN 7191 surf_usm_v(l)%surface_types(m) = pedestrian_category !< default category for wall surface in pedestrian zone 7192 surf_usm_v(l)%albedo(:,m) = -1.0_wp 7193 surf_usm_v(l)%thickness_wall(m) = -1.0_wp 7194 surf_usm_v(l)%thickness_window(m) = -1.0_wp 7195 surf_usm_v(l)%thickness_green(m) = -1.0_wp 7196 surf_usm_v(l)%transmissivity(m) = -1.0_wp 7249 surf_usm_v(l)%surface_types(m) = pedestrian_category !< default category for wall surface in 7250 !<pedestrian zone 7251 surf_usm_v(l)%albedo(:,m) = -1.0_wp 7252 surf_usm_v(l)%thickness_wall(m) = -1.0_wp 7253 surf_usm_v(l)%thickness_window(m) = -1.0_wp 7254 surf_usm_v(l)%thickness_green(m) = -1.0_wp 7255 surf_usm_v(l)%transmissivity(m) = -1.0_wp 7197 7256 ELSE 7198 surf_usm_v(l)%surface_types(m) = usm_par(ii+1,jw,iw)7199 surf_usm_v(l)%albedo(:,m) = usm_val(ij,jw,iw)7200 surf_usm_v(l)%thickness_wall(m) = usm_val(ij+1,jw,iw)7201 surf_usm_v(l)%thickness_window(m) 7202 surf_usm_v(l)%thickness_green(m) 7203 surf_usm_v(l)%transmissivity(m) = 0.0_wp7257 surf_usm_v(l)%surface_types(m) = usm_par(ii+1,jw,iw) 7258 surf_usm_v(l)%albedo(:,m) = usm_val(ij,jw,iw) 7259 surf_usm_v(l)%thickness_wall(m) = usm_val(ij+1,jw,iw) 7260 surf_usm_v(l)%thickness_window(m) = usm_val(ij+1,jw,iw) 7261 surf_usm_v(l)%thickness_green(m) = usm_val(ij+1,jw,iw) 7262 surf_usm_v(l)%transmissivity(m) = 0.0_wp 7204 7263 ENDIF 7205 7264 ELSE IF ( kw <= usm_par(ii+2,jw,iw) ) THEN 7265 ! 7206 7266 !-- wall zone 7207 7267 IF ( usm_par(ii+3,jw,iw) == 0 ) THEN 7208 surf_usm_v(l)%surface_types(m) = wall_category !< default category for wall surface7209 surf_usm_v(l)%albedo(:,m) = -1.0_wp7210 surf_usm_v(l)%thickness_wall(m) = -1.0_wp7211 surf_usm_v(l)%thickness_window(m) 7212 surf_usm_v(l)%thickness_green(m) 7213 surf_usm_v(l)%transmissivity(m) = -1.0_wp7268 surf_usm_v(l)%surface_types(m) = wall_category !< default category for wall surface 7269 surf_usm_v(l)%albedo(:,m) = -1.0_wp 7270 surf_usm_v(l)%thickness_wall(m) = -1.0_wp 7271 surf_usm_v(l)%thickness_window(m) = -1.0_wp 7272 surf_usm_v(l)%thickness_green(m) = -1.0_wp 7273 surf_usm_v(l)%transmissivity(m) = -1.0_wp 7214 7274 ELSE 7215 surf_usm_v(l)%surface_types(m) = usm_par(ii+3,jw,iw)7216 surf_usm_v(l)%albedo(:,m) = usm_val(ij+2,jw,iw)7217 surf_usm_v(l)%thickness_wall(m) = usm_val(ij+3,jw,iw)7218 surf_usm_v(l)%thickness_window(m) 7219 surf_usm_v(l)%thickness_green(m) 7220 surf_usm_v(l)%transmissivity(m) = 0.0_wp7275 surf_usm_v(l)%surface_types(m) = usm_par(ii+3,jw,iw) 7276 surf_usm_v(l)%albedo(:,m) = usm_val(ij+2,jw,iw) 7277 surf_usm_v(l)%thickness_wall(m) = usm_val(ij+3,jw,iw) 7278 surf_usm_v(l)%thickness_window(m) = usm_val(ij+3,jw,iw) 7279 surf_usm_v(l)%thickness_green(m) = usm_val(ij+3,jw,iw) 7280 surf_usm_v(l)%transmissivity(m) = 0.0_wp 7221 7281 ENDIF 7222 7282 ELSE IF ( kw <= usm_par(ii+4,jw,iw) ) THEN 7283 ! 7223 7284 !-- roof zone 7224 7285 IF ( usm_par(ii+5,jw,iw) == 0 ) THEN 7225 surf_usm_v(l)%surface_types(m) = roof_category !< default category for roof surface7226 surf_usm_v(l)%albedo(:,m) = -1.0_wp7227 surf_usm_v(l)%thickness_wall(m) = -1.0_wp7228 surf_usm_v(l)%thickness_window(m) 7229 surf_usm_v(l)%thickness_green(m) 7230 surf_usm_v(l)%transmissivity(m) = -1.0_wp7286 surf_usm_v(l)%surface_types(m) = roof_category !< default category for roof surface 7287 surf_usm_v(l)%albedo(:,m) = -1.0_wp 7288 surf_usm_v(l)%thickness_wall(m) = -1.0_wp 7289 surf_usm_v(l)%thickness_window(m) = -1.0_wp 7290 surf_usm_v(l)%thickness_green(m) = -1.0_wp 7291 surf_usm_v(l)%transmissivity(m) = -1.0_wp 7231 7292 ELSE 7232 surf_usm_v(l)%surface_types(m) = usm_par(ii+5,jw,iw)7233 surf_usm_v(l)%albedo(:,m) = usm_val(ij+4,jw,iw)7234 surf_usm_v(l)%thickness_wall(m) = usm_val(ij+5,jw,iw)7235 surf_usm_v(l)%thickness_window(m) 7236 surf_usm_v(l)%thickness_green(m) 7237 surf_usm_v(l)%transmissivity(m) = 0.0_wp7293 surf_usm_v(l)%surface_types(m) = usm_par(ii+5,jw,iw) 7294 surf_usm_v(l)%albedo(:,m) = usm_val(ij+4,jw,iw) 7295 surf_usm_v(l)%thickness_wall(m) = usm_val(ij+5,jw,iw) 7296 surf_usm_v(l)%thickness_window(m) = usm_val(ij+5,jw,iw) 7297 surf_usm_v(l)%thickness_green(m) = usm_val(ij+5,jw,iw) 7298 surf_usm_v(l)%transmissivity(m) = 0.0_wp 7238 7299 ENDIF 7239 7300 ELSE 7240 !7241 7301 WRITE(9,*) 'Problem reading USM data:' 7242 7302 WRITE(9,*) l,i,j,kw,get_topography_top_index_ji( j, i, 's' ) … … 7247 7307 WRITE(9,*) kw,roof_height_limit,wall_category,roof_category 7248 7308 FLUSH(9) 7309 ! 7249 7310 !-- supply the default category 7250 7311 IF ( kw <= roof_height_limit ) THEN … … 7253 7314 surf_usm_v(l)%surface_types(m) = roof_category !< default category for wall surface in roof zone 7254 7315 END IF 7255 surf_usm_v(l)%albedo(:,m) = -1.0_wp7256 surf_usm_v(l)%thickness_wall(m) = -1.0_wp7316 surf_usm_v(l)%albedo(:,m) = -1.0_wp 7317 surf_usm_v(l)%thickness_wall(m) = -1.0_wp 7257 7318 surf_usm_v(l)%thickness_window(m) = -1.0_wp 7258 surf_usm_v(l)%thickness_green(m) = -1.0_wp7259 surf_usm_v(l)%transmissivity(m) = -1.0_wp7319 surf_usm_v(l)%thickness_green(m) = -1.0_wp 7320 surf_usm_v(l)%transmissivity(m) = -1.0_wp 7260 7321 ENDIF 7261 7322 ! … … 7270 7331 ENDDO 7271 7332 IF ( ip == -99999 ) THEN 7333 ! 7272 7334 !-- wall category not found 7273 7335 WRITE (9, "(A,I7,A,3I5)") 'wall category ', it, & … … 7282 7344 ENDDO 7283 7345 IF ( ip == -99999 ) THEN 7346 ! 7284 7347 !-- default wall category not found 7285 7348 WRITE (9, "(A,I5,A,3I5)") 'Default wall category', category, ' not found!' … … 7405 7468 ENDDO 7406 7469 ENDIF 7407 7470 ! 7408 7471 !-- coordinate not found 7409 7472 isurfl = -1 … … 7423 7486 SUBROUTINE usm_read_wall_temperature 7424 7487 7425 INTEGER(iwp) :: i, j, k, d, ii, iline 7488 INTEGER(iwp) :: i, j, k, d, ii, iline !> running indices 7426 7489 INTEGER(iwp) :: isurfl 7427 7490 REAL(wp) :: rtsurf … … 7429 7492 7430 7493 7431 7432 7433 7494 DO ii = 0, io_blocks-1 7434 7495 IF ( ii == io_group ) THEN 7435 7496 ! 7436 7497 !-- open wall temperature file 7437 7498 OPEN( 152, file='WALL_TEMPERATURE'//coupling_char, action='read', & … … 7454 7515 CALL message( 'usm_read_wall_temperature', 'PA0521', 1, 2, 0, 6, 0 ) 7455 7516 ENDIF 7456 7517 ! 7457 7518 !-- assign temperatures 7458 7519 IF ( d == 0 ) THEN … … 7498 7559 !> with many simplifications and adjustments. 7499 7560 !> TODO better description 7561 !> No calculation of window surface temperatures during spinup to increase 7562 !> maximum possible timstep 7500 7563 !------------------------------------------------------------------------------! 7501 7564 SUBROUTINE usm_surface_energy_balance( spinup ) … … 7504 7567 IMPLICIT NONE 7505 7568 7506 INTEGER(iwp) :: i, j, k, l, m !< running indices7569 INTEGER(iwp) :: i, j, k, l, m !< running indices 7507 7570 7508 7571 INTEGER(iwp) :: i_off !< offset to determine index of surface element, seen from atmospheric grid point, for x … … 7512 7575 LOGICAL :: spinup !true during spinup 7513 7576 7514 REAL(wp) :: stend_wall 7577 REAL(wp) :: stend_wall !< surface tendency 7515 7578 7516 7579 REAL(wp) :: stend_window !< surface tendency … … 7526 7589 REAL(wp) :: f_shf_window !< factor for shf_eb window 7527 7590 REAL(wp) :: f_shf_green !< factor for shf_eb green wall 7528 REAL(wp) :: lambda_surface !< current value of lambda_surface (heat conductivity between air and wall) 7529 REAL(wp) :: lambda_surface_window !< current value of lambda_surface (heat conductivity between air and window) 7530 REAL(wp) :: lambda_surface_green !< current value of lambda_surface (heat conductivity between air and greeb wall) 7591 REAL(wp) :: lambda_surface !< current value of lambda_surface (heat conductivity 7592 !<between air and wall) 7593 REAL(wp) :: lambda_surface_window !< current value of lambda_surface (heat conductivity 7594 !< between air and window) 7595 REAL(wp) :: lambda_surface_green !< current value of lambda_surface (heat conductivity 7596 !< between air and greeb wall) 7531 7597 7532 7598 REAL(wp) :: dtime !< simulated time of day (in UTC) 7533 7599 INTEGER(iwp) :: dhour !< simulated hour of day (in UTC) 7534 7600 REAL(wp) :: acoef !< actual coefficient of diurnal profile of anthropogenic heat 7535 REAL(wp) :: f1, & !< resistance correction term 17536 f2, & !< resistance correction term 27537 f3, & !< resistance correction term 37538 e, & !< water vapour pressure7539 e_s, & !< water vapour saturation pressure7540 e_s_dt, & !< derivate of e_s with respect to T7541 tend, & !< tendency7542 dq_s_dt, & !< derivate of q_s with respect to T7543 f_qsws, & !< factor for qsws7544 f_qsws_veg, & !< factor for qsws_veg7545 f_qsws_liq, & !< factor for qsws_liq7546 m_liq_max, & !< maxmimum value of the liq. water reservoir7547 qv1, & !< specific humidity at first grid level7548 m_max_depth = 0.0002_wp, & !Maximum capacity of the water reservoir (m)7549 rho_lv, & 7550 drho_l_lv, & 7551 q_s 7601 REAL(wp) :: f1, & !< resistance correction term 1 7602 f2, & !< resistance correction term 2 7603 f3, & !< resistance correction term 3 7604 e, & !< water vapour pressure 7605 e_s, & !< water vapour saturation pressure 7606 e_s_dt, & !< derivate of e_s with respect to T 7607 tend, & !< tendency 7608 dq_s_dt, & !< derivate of q_s with respect to T 7609 f_qsws, & !< factor for qsws 7610 f_qsws_veg, & !< factor for qsws_veg 7611 f_qsws_liq, & !< factor for qsws_liq 7612 m_liq_max, & !< maxmimum value of the liq. water reservoir 7613 qv1, & !< specific humidity at first grid level 7614 m_max_depth = 0.0002_wp, & !< Maximum capacity of the water reservoir (m) 7615 rho_lv, & !< frequently used parameter for green layers 7616 drho_l_lv, & !< frequently used parameter for green layers 7617 q_s !< saturation specific humidity 7552 7618 7553 7619 ! … … 7602 7668 !-- pt, us, ts are not available for the prognostic time step, 7603 7669 !-- data from the last time step is used here. 7604 7670 ! 7605 7671 !-- Workaround: use single r_a as stability is only treated for the 7606 7672 !-- average temperature … … 7621 7687 surf_usm_h%r_a(m) = 1.0_wp 7622 7688 IF ( surf_usm_h%r_a_green(m) < 1.0_wp ) & 7623 surf_usm_h%r_a_green(m) = 1.0_wp7689 surf_usm_h%r_a_green(m) = 1.0_wp 7624 7690 IF ( surf_usm_h%r_a_window(m) < 1.0_wp ) & 7625 7691 surf_usm_h%r_a_window(m) = 1.0_wp … … 7631 7697 surf_usm_h%r_a(m) = 300.0_wp 7632 7698 IF ( surf_usm_h%r_a_green(m) > 300.0_wp ) & 7633 surf_usm_h%r_a_green(m) = 300.0_wp7699 surf_usm_h%r_a_green(m) = 300.0_wp 7634 7700 IF ( surf_usm_h%r_a_window(m) > 300.0_wp ) & 7635 7701 surf_usm_h%r_a_window(m) = 300.0_wp 7636 7702 7637 7703 ! 7638 7704 !-- factor for shf_eb 7639 7705 f_shf = rho_cp / surf_usm_h%r_a(m) … … 7641 7707 f_shf_green = rho_cp / surf_usm_h%r_a_green(m) 7642 7708 7643 !*************************************************************************************** 7644 if (surf_usm_h%frac(ind_pav_green,m).gt.0.0_wp) then 7645 !-- Adapted from LSM:7646 !-- Second step: calculate canopy resistance r_canopy7647 !-- f1-f3 here are defined as 1/f1-f3 as in ECMWF documentation7709 7710 IF ( surf_usm_h%frac(ind_pav_green,m) > 0.0_wp ) THEN 7711 !-- Adapted from LSM: 7712 !-- Second step: calculate canopy resistance r_canopy 7713 !-- f1-f3 here are defined as 1/f1-f3 as in ECMWF documentation 7648 7714 7649 !-- f1: correction for incoming shortwave radiation (stomata close at7650 !-- night)7651 f1 = MIN( 1.0_wp, ( 0.004_wp * surf_usm_h%rad_sw_in(m) + 0.05_wp ) / &7652 (0.81_wp * (0.004_wp * surf_usm_h%rad_sw_in(m) &7653 + 1.0_wp)) )7654 ! 7655 !-- f2: correction for soil moisture availability to plants (the7656 !-- integrated soil moisture must thus be considered here)7657 !-- f2 = 0 for very dry soils7658 m_total = 0.0_wp7659 DO k = nzb_wall, nzt_wall+17660 m_total = m_total + rootfr_h(nzb_wall,m) &7661 * MAX(swc_h(nzb_wall,m),wilt_h(nzb_wall,m))7662 ENDDO7663 7664 IF ( m_total > wilt_h(nzb_wall,m) .AND. m_total < fc_h(nzb_wall,m) ) THEN7665 f2 = ( m_total - wilt_h(nzb_wall,m) ) / (fc_h(nzb_wall,m) - wilt_h(nzb_wall,m) )7666 ELSEIF ( m_total >= fc_h(nzb_wall,m) ) THEN7667 f2 = 1.0_wp7668 ELSE7669 f2 = 1.0E-20_wp7670 ENDIF7671 7672 ! 7673 !-- Calculate water vapour pressure at saturation7674 e_s = 0.01_wp * 610.78_wp * EXP( 17.269_wp * ( t_surf_green_h(m) &7675 - 273.16_wp ) / ( t_surf_green_h(m) - 35.86_wp ) )7676 ! 7677 !-- f3: correction for vapour pressure deficit7678 IF ( surf_usm_h%g_d(m) /= 0.0_wp ) THEN7679 ! 7680 !-- Calculate vapour pressure7681 e = qv1 * surface_pressure / ( qv1 + 0.622_wp )7682 f3 = EXP ( - surf_usm_h%g_d(m) * (e_s - e) )7683 ELSE7684 f3 = 1.0_wp7685 ENDIF7686 7687 ! 7688 !-- Calculate canopy resistance. In case that c_veg is 0 (bare soils),7689 !-- this calculation is obsolete, as r_canopy is not used below.7690 !-- To do: check for very dry soil -> r_canopy goes to infinity7691 surf_usm_h%r_canopy(m) = surf_usm_h%r_canopy_min(m) /&7715 !-- f1: correction for incoming shortwave radiation (stomata close at 7716 !-- night) 7717 f1 = MIN( 1.0_wp, ( 0.004_wp * surf_usm_h%rad_sw_in(m) + 0.05_wp ) / & 7718 (0.81_wp * (0.004_wp * surf_usm_h%rad_sw_in(m) & 7719 + 1.0_wp)) ) 7720 ! 7721 !-- f2: correction for soil moisture availability to plants (the 7722 !-- integrated soil moisture must thus be considered here) 7723 !-- f2 = 0 for very dry soils 7724 m_total = 0.0_wp 7725 DO k = nzb_wall, nzt_wall+1 7726 m_total = m_total + rootfr_h(nzb_wall,m) & 7727 * MAX(swc_h(nzb_wall,m),wilt_h(nzb_wall,m)) 7728 ENDDO 7729 7730 IF ( m_total > wilt_h(nzb_wall,m) .AND. m_total < fc_h(nzb_wall,m) ) THEN 7731 f2 = ( m_total - wilt_h(nzb_wall,m) ) / (fc_h(nzb_wall,m) - wilt_h(nzb_wall,m) ) 7732 ELSEIF ( m_total >= fc_h(nzb_wall,m) ) THEN 7733 f2 = 1.0_wp 7734 ELSE 7735 f2 = 1.0E-20_wp 7736 ENDIF 7737 7738 ! 7739 !-- Calculate water vapour pressure at saturation 7740 e_s = 0.01_wp * 610.78_wp * EXP( 17.269_wp * ( t_surf_green_h(m) & 7741 - 273.16_wp ) / ( t_surf_green_h(m) - 35.86_wp ) ) 7742 ! 7743 !-- f3: correction for vapour pressure deficit 7744 IF ( surf_usm_h%g_d(m) /= 0.0_wp ) THEN 7745 ! 7746 !-- Calculate vapour pressure 7747 e = qv1 * surface_pressure / ( qv1 + 0.622_wp ) 7748 f3 = EXP ( - surf_usm_h%g_d(m) * (e_s - e) ) 7749 ELSE 7750 f3 = 1.0_wp 7751 ENDIF 7752 7753 ! 7754 !-- Calculate canopy resistance. In case that c_veg is 0 (bare soils), 7755 !-- this calculation is obsolete, as r_canopy is not used below. 7756 !-- To do: check for very dry soil -> r_canopy goes to infinity 7757 surf_usm_h%r_canopy(m) = surf_usm_h%r_canopy_min(m) / & 7692 7758 ( surf_usm_h%lai(m) * f1 * f2 * f3 + 1.0E-20_wp ) 7693 7759 7694 7760 !