Changeset 183 for palm


Ignore:
Timestamp:
Aug 4, 2008 3:39:12 PM (16 years ago)
Author:
letzel
Message:
  • mrun: (optionally) link local input
  • diffusion_s: bugfix for calculation of horizontal fluxes at vertical walls
  • NCL scripts in trunk/SCRIPTS/NCL updated
Location:
palm/trunk
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/SCRIPTS/NCL/spectra.ncl

    r178 r183  
    1313      delete(parameter@_FillValue)
    1414   else
    15       if (isfilepresent("~/palm/current_version/trunk/SRIPTS/NCL/.ncl_preferences")) then
    16          parameter = asciiread("~/palm/current_version/trunk/SRIPTS/NCL/.ncl_preferences",129,"string")
     15      if (isfilepresent("~/palm/current_version/trunk/SCRIPTS/NCL/.ncl_preferences")) then
     16         parameter = asciiread("~/palm/current_version/trunk/SCRIPTS/NCL/.ncl_preferences",129,"string")
    1717         delete(parameter@_FillValue)
    1818      else
    1919         print(" ")
    20          print("'.ncl_preferences' does not exist in '~/palm/current_version/trunk/SRIPTS/NCL/'")
     20         print("'.ncl_preferences' does not exist in '~/palm/current_version/trunk/SCRIPTS/NCL/'")
    2121         print(" ")
    2222         exit
     
    372372   if (sort .EQ. "time")
    373373      plot  = new(dim*dimz,graphic)
    374       plot_ = new(dim*dimz,graphic)
    375374      np=dimt
    376375      res@lgTitleString = "Time [h]"
    377376   else
    378377      plot  = new(dim*dimt,graphic)
    379       plot_ = new(dim*dimt,graphic)
    380378      np=dimz
    381379      res@lgTitleString = "Height [m]" 
     
    455453               res@trXMaxF = max(x_axis)
    456454               res@gsnLeftString      = vNam(varn)
    457                res@gsnRightString     = "Height: "+z(p)+"m"
     455               res@gsnRightString     = "Height = "+z(p)+"m"
    458456               if (norm .NE. 1)then
    459457                  res@tiYAxisString      = vNam(varn)+" / "+norm
     
    482480                  res@trXMaxF = max(x_axis)
    483481                  res@gsnLeftString      = vNam(varn)
    484                   res@gsnRightString     = "Time: "+legend_label(p)+"h"
     482                  res@gsnRightString     = "Time = "+legend_label(p)+"h"
    485483                  if (norm .NE. 1)then
    486484                     res@tiYAxisString      = vNam(varn)+" / "+norm
     
    521519   resP@txFontHeightF          = 0.014
    522520
    523    do m=0,n-1
    524       plot_(m)=plot(n-1-m)
    525    end do
    526 
    527521   if (format_out .EQ. "eps" .OR. format_out .EQ. "epsi") then
    528       gsn_panel(wks,plot_,(/n,1/),resP)
     522      gsn_panel(wks,plot,(/n,1/),resP)
    529523   else   
    530524      do i = 0,n-1, no_lines*no_columns
    531525         if( (i+no_lines*no_columns) .gt. (n-1)) then
    532             gsn_panel(wks,plot_(i:n-1),(/no_lines,no_columns/),resP)
     526            gsn_panel(wks,plot(i:n-1),(/no_lines,no_columns/),resP)
    533527         else
    534            gsn_panel(wks,plot_(i:i+no_lines*no_columns-1),(/no_lines,no_columns/),resP)
     528           gsn_panel(wks,plot(i:i+no_lines*no_columns-1),(/no_lines,no_columns/),resP)
    535529         end if
    536530      end do
  • palm/trunk/SCRIPTS/mrun

    r182 r183  
    147147     #                     adjustments for lcxt4 (Bergen Center for Computational
    148148     #                     Science)
     149     # 22/05/08 - Marcus - If environment variable link_local_input is set to
     150     #                     true, mrun tries "ln -f" on local input and resorts
     151     #                     to "cp" or "cp -r" on error
    149152     # 27/05/08 - Siggi  - PATH is set to PALM_BIN everywhere (missing so far)
    150153     # 14/07/08 - Siggi  - adjustments for lcsgih
    151 
    152154 
    153155    # VARIABLENVEREINBARUNGEN + DEFAULTWERTE
     
    196198 job_on_file=""
    197199 keep_data_from_previous_run=false
     200 link_local_input=false
    198201 link_local_output=false
    199202 localhost_realname=$(hostname)
     
    25392542                printf "\n      providing $numprocs files for the respective processors"
    25402543                mkdir  ${localin[$i]}
    2541                 cp -r  ${absnamein[$i]}/*  ${localin[$i]}
     2544                if [[ $link_local_input = true ]]
     2545                    then
     2546                    printf "      using ln -f\n"
     2547                    cd ${absnamein[$i]}
     2548                    for file in $(ls *)
     2549                      do
     2550                      ln -f $file  ${localin[$i]}
     2551                    done
     2552                    cd $TEMPDIR
     2553                fi
     2554                # If "ln -f" fails of if "$link_local_input = false" do a normal "cp -r"
     2555                if [[ ! -f "${localin[$i]}/_0000" ]]
     2556                    then
     2557                    if [[ $link_local_input = true ]]
     2558                        then
     2559                        printf "      ln failed for .../_0000, using cp...\n"
     2560                    fi
     2561                    cp -r  ${absnamein[$i]}/*  ${localin[$i]}
     2562                fi
    25422563
    25432564             else
    25442565                   # BEREITSTELLUNG AUF EINPROZESSORRECHNERN
    2545                 cp  ${absnamein[$i]}  ${localin[$i]}
     2566                if [[ $link_local_input = true ]]
     2567                then
     2568                    printf "      using ln -f\n"
     2569                    ln -f  ${absnamein[$i]}  ${localin[$i]}
     2570                fi
     2571                # If "ln -f" fails of if "$link_local_input = false" do a normal "cp"
     2572                if [[ ! -f "${localin[$i]}" ]]
     2573                then
     2574                    if [[ $link_local_input = true ]]
     2575                    then
     2576                        printf "      ln failed, using cp...\n"
     2577                    fi
     2578                    cp  ${absnamein[$i]}  ${localin[$i]}
     2579                fi
    25462580             fi
    25472581          fi
  • palm/trunk/SOURCE/CURRENT_MODIFICATIONS

    r181 r183  
    8181Errors:
    8282------
     83Bugfix: calculation of horizontal fluxes at vertical walls (diffusion_s)
     84
    8385Bugfix: zero assignments to tendency arrays in case of restarts (init_3d_model)
    8486
     
    117119
    118120
    119 check_parameters, flow_statistics, init_dvrp, init_3d_model, local_stop, plant_canopy_model, poismg, pres, read_3d_binary, user_interface, write_3d_binary
     121check_parameters, diffusion_s, flow_statistics, init_dvrp, init_3d_model, local_stop, plant_canopy_model, poismg, pres, read_3d_binary, user_interface, write_3d_binary
    120122
  • palm/trunk/SOURCE/diffusion_s.f90

    r139 r183  
    44! Actual revisions:
    55! -----------------
    6 !
     6! bugfix: calculation of fluxes at vertical surfaces
    77!
    88! Former revisions:
     
    9090
    9191                   tend(k,j,i) = tend(k,j,i)                                  &
    92                                           + 0.5 * ( fwxp(j,i) *               &
     92                                                + ( fwxp(j,i) * 0.5 *         &
    9393                        ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
    9494                        + ( 1.0 - fwxp(j,i) ) * wall_s_flux(1)                &
    95                                                    -fwxm(j,i) *               &
     95                                                   -fwxm(j,i) * 0.5 *         &
    9696                        ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
    9797                        + ( 1.0 - fwxm(j,i) ) * wall_s_flux(2)                &
    9898                                                  ) * ddx2                    &
    99                                           + 0.5 * ( fwyp(j,i) *               &
     99                                                + ( fwyp(j,i) * 0.5 *         &
    100100                        ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
    101101                        + ( 1.0 - fwyp(j,i) ) * wall_s_flux(3)                &
    102                                                    -fwym(j,i) *               &
     102                                                   -fwym(j,i) * 0.5 *         &
    103103                        ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
    104104                        + ( 1.0 - fwym(j,i) ) * wall_s_flux(4)                &
     
    200200
    201201             tend(k,j,i) = tend(k,j,i)                                        &
    202                                           + 0.5 * ( fwxp(j,i) *               &
     202                                                + ( fwxp(j,i) * 0.5 *         &
    203203                        ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
    204204                        + ( 1.0 - fwxp(j,i) ) * wall_s_flux(1)                &
    205                                                    -fwxm(j,i) *               &
     205                                                   -fwxm(j,i) * 0.5 *         &
    206206                        ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
    207207                        + ( 1.0 - fwxm(j,i) ) * wall_s_flux(2)                &
    208208                                                  ) * ddx2                    &
    209                                           + 0.5 * ( fwyp(j,i) *               &
     209                                                + ( fwyp(j,i) * 0.5 *         &
    210210                        ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
    211211                        + ( 1.0 - fwyp(j,i) ) * wall_s_flux(3)                &
    212                                                    -fwym(j,i) *               &
     212                                                   -fwym(j,i) * 0.5 *         &
    213213                        ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
    214214                        + ( 1.0 - fwym(j,i) ) * wall_s_flux(4)                &
Note: See TracChangeset for help on using the changeset viewer.