damyarou

python, GMT などのプログラム

macOS Big Sur インストール

2020年11月13日 MacBook Pro (Retina, 13-inch, Mid 2014)

今日11月13日、気がついたらMacbookにアップデートの知らせが。なんだろうと思ったらBig Surのアップデートの通知でした。新しもの好きな私は、早速アップデートを実行。流石にメインで使っているiMacはCatalinaのままいじらず。アップデートにかかった時間は正確にはわかりませんが、ダウンロードから再起動まで50〜60分程度だと思います。特にトラブルもなくスムースにアップデート完了。

Norton OK、ネットOK、4KディスプレイOK、VS code OK、Jupyter notebook OKでとりあえず普段使っているものは動きます（詳細まではまだ確認していません）。確かにUIは少し丸くなり、かわいくなった感じです。

f:id:damyarou:20201122081242j:plain

ファインダーの文字なども見やすくなった感じ。まだ15分ほどいじっただけなので、挙動についての報告はありません。「とりあえず早いものがち」という感じの記事でした。

なお、古めのMacbook ProをBig Surにアップデートすると修復不能のトラブルが発生するとの情報もあります。

古めのMacBook Pro、Big Surアップデートで文鎮化したとの報告多数

アップル、Big Surアップデートで文鎮化MacBook Proの救出方法を説明

しかしながら、私の場合は問題なくインストールでき、動いていますので、運が良かった？

2020年11月14日 iMac (Retina 4K, 21.5-inch, 2017)

iMac (Retina 4K, 21.5-inch, 2017)にもBig Surいれました。 インストールはダウンロード (12.18GB) 含めて45分弱。インストール中のトラブルは無し、現時点では不具合なく使えています。

f:id:damyarou:20201122081326j:plain

以　上

iPad Air 4 のためのグッズ

iPad Air 4のためのグッズが少し揃ってきたので紹介します。こうして並べてみると、結構お金かけているなあ。iPad Airでやりたいことは、主に手書きメモとpdf文書のチェック（赤書き・マーキング）なので、今の所キーボードの購入予定は無しです。キーボードを使うことはMacで。

（新型Macbook AirとMacbook Proの予約が開始になったようです。買い換えるか考案中。もう少しまとうかなあ）

  • iPad Air 4: Apple Store 2020/11/01 注文、2020/11/08到着 (¥62,800)
  • Apple Pencil（第２世代）: Apple Store 2020/11/01 注文、2020/11/02到着 (¥14,500)
  • ペーパーライクフィルタ Bellemond ペン先の消耗を抑える/ケント紙: Amazon 2020/11/04注文、2020/11/06到着 (¥1,698)
  • ESR iPad Air 4 ケース2020 iPad 10.9インチ（ネイビーブルー）: Amazon 2020/11/12注文、2020/11/13到着 (¥2,099)
  • GoodNotes 5: Apple Store 2020/11/10購入 (¥980)

以下にESRのケースの写真を紹介します。色はネイビーブルー。iPadがスカイブルーなので、青系にしてみました。カメラの穴の位置・形ともOKです。Apple Pencilも充電しながら収納できます。Tatch IDも先に登録しておけば、ケースをつけた後は特に使いづらいとは感じません。

f:id:damyarou:20201122080856j:plain

f:id:damyarou:20201122080930j:plain

f:id:damyarou:20201122080956j:plain

メモアプリとしてGoodNotes 5を購入。ノートアプリも色々あり、どれにするか迷うところですが、ネット上で比較的評判がよく、やりたいこと（手書きメモ作成とPDFへの赤入れ作業程度）を、便利にこなせそうということで、GoodNotes 5にしました。下の写真はSplit viewでGodNotesを開いているところ。左側に取り込んだpdfを見ながら右でメモを取れます。pdfにもコメントを入れたり蛍光ペンマークできます。pdf書類そのものはiMacあるいはMacbookで作りますが、これをiCloud経由もしくはAirDropでiPadに取り込みます。iPad AirはUSB Type-C端子を持っているので、外付け記憶デバイスも使ってみたいのですが、もう少し勉強してからにします。iPad Airで２つ資料を並べると少し画面が小さい気もしますが、作っているのはメモなので、「こういうもの」として納得しています。10.9-inchという大きさはA4用紙より一回り小さいB5用紙程度であり、他の書類や書籍と一緒に持ち歩くには、この程度の大きさがちょうどいいと思っています。

f:id:damyarou:20201122081034j:plain

こういうものは、とのかく触って、やりたいことがデキるかを自分で確認していくことが使いこなすための秘訣なのでしょう。これからも色々使っていこうと思います。

以　上

iPad Air 4 購入

iPad Air 4 注文

2020年11月1日、iPad Air 10.9インチを、Apple pencil とともに注文しました。iPad Air は、ストレージ64GB、スカイブルーのWi-Fiモデルです。Amazon か Apple store か、どちらに注文するか迷いましたが、Amazon では入荷が不安定のようなので、本家である Apple store を利用しました。Apple pencil は在庫があるということで、翌日11月2日に手元に届きました。iPad 本体は当初は11月20日前後のお届けということでしたが、その後11月8日配達ということでメールが来ました。Apple pencil はまだ開封せず、iPad 本体と主に開封・写真撮影を行うつもりです。

今使っているタブレットは iPad mini 4 です。iPad mini を購入する前は、以前より持っていた Macbook Pro (Retina, 13-inch, Mid 2014) を持ち歩き、新幹線の中などで Web 閲覧を行ったりしていたのですが、Macbook Pro を開きながら Web 閲覧というのも大袈裟すぎるし恥ずかしいので、iPad mini を購入し、新幹線の中などで使うようにしました。ちなみに私は新幹線通勤をしており、通勤時間は片道２時間、そのうち１時間は新幹線の中です。現在は出社は週１回にし、その他は在宅勤務としています。iPad mini は新幹線の中の他、家にいるときでも寝る前や朝起きがけにニュースを見たり音楽を聞いたりメールチェックなどに使っています。パソコンと違い、すぐに使えるところがとても便利です。

iPad mini 4 もそれなりに活躍していたのですが、なぜ新しく出た iPad Air を買うことにしたのでしょうか？主に着目した点は２つあります。１つめは、Apple pencil 第二世代に対応したこと、２つめは USB-C を積んだことです。iPad mini 4 もたまーにですが、ソフトキーボードでメモしたりという使い方もしていましたが、ネットでも手書きメモの話題がよく目に付き、自分でもApple pencil を使って手書きメモやスケッチをやってみたいと思うようになりました。これまで、メモはいらなくなったA4印刷物を半分に切って裏面を使っていましたが、空白部がある程度以下になると捨てていたのですが、あとから「こないだメモしたあれはどうしたっけ？？？」ということも結構あり、メモをスマートに保存しておく必要性も感じています。また USB-C に対応したため、外付けストレージ経由でiMacとの連携で運用できれば使いみちが広がるかもしれないと思っているところです。なお、新しいiPad Airにイヤホンジャックが無いのが痛いですが、そこはUSC端子を使ってなんとかできそうです。ちなみに Macbook Pro Mid 2014 は USB 端子は A のみですが、そこはアダプタを使えばこちらとも連携できそう（やってみなければわかりませんが）。ただし、Macbook pro は次に最新型が出たら買い換えようかとも思っています。

iPad Air 4 到着

2020年11月8日、昼前に予定通り11月1日に注文したiPad Air 4が届きました。下の写真右がiPad Air、真ん中が先に届いていたApple pencil、左がこれも先に注文していたペーパーライクフィルム（保護フィルム）。BELLEMONDという会社の「ペン先の消耗を抑える／ケント紙」というやつです。保護フィルムはどれにするか迷いましたが結局は直感で決めました。

f:id:damyarou:20201122080142j:plain

次がiPhone 8、iPad mini 4、iPad Air 4の集合写真です。

f:id:damyarou:20201122080313j:plain

iPad Airは保護フィルムを貼ってある状態です。画面の明るさはデフォルトのままです。透明つるつるの保護フィルムを貼ってあるiPad mini 4と比べると少し暗い感じかな。表面もザラザラ感があります。Apple pencilも試し書きしてみましたが、まだ慣れてなくて硬いものの上に書いている感じ。いずれ慣れるでしょう。

f:id:damyarou:20201122080343j:plain

上はアップルマーク側の写真です。iPhone 8は白色、iPad mini 4は金色、iPad Air 4は青色です。iPad Air は、下に示す銀色のMacbook Pro 13-inchと比べると明らかに色が違う。青と言われれば青ですよね（私は目が悪いため色の識別には苦労します）。

f:id:damyarou:20201122080413j:plain

こうやって機械を並べて写真を撮っているとそれだけで楽しくなってきます。「デキる男」になったような錯覚に陥ります。が、、、iPad Airを買ったのはメモ（ある計算をするためのメモを要らなくなった紙の裏に鉛筆描きしたもの）をiPad Air上に描ききちんととっておくためなので、本来の目的を達成できるよう、使い方を勉強します。

iPad Airのケースはまだ仕入れていないため、キズをつけないよう早めに専用ケースを手に入れようと思います。これから色々いじるのが楽しみです。

以　上

Stability analysis of large dam non-overflow section

Code is shown below.

# Concrete gravity dam stability analysis
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.ticker as tick


#global variants
wc=23.5 # unit weight of concrete (kN/m3) 
ww=10.0 # unit weight of water (kN/m3) 
ws=12.0 # unit weight of sedimrntation (submerged) (kN/m3) 
ce=0.5  # sedimentation pressure coefficient
kh=0.15 # seismic coefficient
cc=3200 # cohesion of foundation rock (kN/m2) 
ff=1.20 # internal friction coefficient of foundation rock


def barrow(x1,y1,x2,y2):
    col='#333333'
    arst='<->,head_width=3,head_length=10'
    aprop=dict(arrowstyle=arst,connectionstyle='arc3',facecolor=col,edgecolor=col,shrinkA=0,shrinkB=0,lw=1)
    plt.annotate('',
        xy=(x1,y1), xycoords='data',
        xytext=(x2,y2), textcoords='data', fontsize=0, arrowprops=aprop)

    
def drawfig(tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m,fpara,fnameF):
    fsz=12
    xmin=fpara[0]
    xmax=fpara[1]
    dx  =fpara[2]
    ymin=fpara[3]
    ymax=fpara[4]
    dy  =fpara[5]
    iw=12
    ih=iw*(ymax-ymin)/(xmax-xmin)
    plt.figure(figsize=(iw,ih),facecolor='w')
    plt.rcParams['font.size']=fsz
    plt.rcParams['font.family']='sans-serif'

    plt.xlim([xmin,xmax])
    plt.ylim([ymin,ymax])
    plt.xlabel('Distance (m)')
    plt.ylabel('Elevation (EL.m)')
    plt.xticks(np.arange(xmin,xmax+dx,dx))
    plt.yticks(np.arange(ymin,ymax+dy,dy))
    plt.gca().spines['right'].set_visible(False)
    plt.gca().spines['top'].set_visible(False)
    plt.gca().yaxis.set_ticks_position('left')
    plt.gca().xaxis.set_ticks_position('bottom')
    plt.gca().set_aspect('equal',adjustable='box')
    plt.title(tstr,loc='left',fontsize=fsz)

    plt.fill([xmin+5,xmax-5,xmax-5,xmin+5],[ELB,ELB,ELB-5,ELB-5],fill=False, hatch='///',color='#aaaaaa',lw=0)
    plt.fill([xmin+5,0,0,xmin+5],[UWL,UWL,UWL-10,UWL-10],color='#00ffff',alpha=0.3,lw=0)
    plt.fill([xmin+5,0,0,xmin+5],[USL,USL,USL-10,USL-10],fill=False, hatch='...',color='#aaaaaa',lw=0)
    x1=n*(ELT-(DSL-10)); x2=n*(ELT-DSL); x3=xmax-5; x4=x3
    y1=DSL-10; y2=DSL; y3=y2; y4=y1
    plt.fill([x1,x2,x3,x4],[y1,y2,y3,y4],fill=False, hatch='...',color='#aaaaaa',lw=0)
    x1=n*(ELT-(DWL-10)); x2=n*(ELT-DWL); x3=xmax-5; x4=x3
    y1=DWL-10; y2=DWL; y3=y2; y4=y1
    plt.fill([x1,x2,x3,x4],[y1,y2,y3,y4],color='#00ffff', alpha=0.3,lw=0)
    
    x=np.array([-m*(ELF-ELB),0,0,tw,tw,n*(ELT-ELB)])
    y=np.array([ELB,ELF,ELT,ELT,ELT-tw/n,ELB])    
    plt.fill(x,y,facecolor='#eeeeee',edgecolor='#000000',lw=2)
    plt.plot([0,tw],[ELT,ELT-tw/n],'--',color='#000000',lw=1)
    
    plt.plot([xmin+5,xmax-5],[ELB,ELB],'-',color='#000000',lw=2)
    plt.plot([xmin+5,0],[UWL,UWL],'-',color='#000000',lw=1)
    plt.plot([xmin+5,0],[USL,USL],'-',color='#000000',lw=1)
    plt.plot([n*(ELT-DWL),xmax-5],[DWL,DWL],'-',color='#000000',lw=1)
    plt.plot([n*(ELT-DSL),xmax-5],[DSL,DSL],'-',color='#000000',lw=1)
    plt.text(xmin+5,ELB,'EL.{0:.3f}'.format(ELB),ha='left',va='bottom',fontsize=fsz)
    plt.text(xmax-5,ELB,'EL.{0:.3f}'.format(ELB),ha='right',va='bottom',fontsize=fsz)
    plt.text(xmin+5,UWL,'UWL:EL.{0:.3f}'.format(UWL),ha='left',va='bottom',fontsize=fsz)
    plt.text(xmin+5,USL,'USL:EL.{0:.3f}'.format(USL),ha='left',va='bottom',fontsize=fsz)
    plt.text(xmax-5,DWL,'DWL:EL.{0:.3f}'.format(DWL),ha='right',va='bottom',fontsize=fsz)
    plt.text(xmax-5,DSL,'DSL:EL.{0:.3f}'.format(DSL),ha='right',va='top',fontsize=fsz)
    plt.plot([tw,tw+30],[ELT,ELT],'-',color='#000000',lw=1)
    plt.text(tw+30,ELT,'EL.{0:.3f}'.format(ELT),ha='right',va='bottom',fontsize=fsz)
    xs=0.5*n*(ELT-ELB); ys=0.5*(ELT+ELB); ang=np.degrees(np.arctan(1/n))
    plt.text(xs+5,ys,'1:{0:.2f}'.format(n),va='center',ha='center',rotation=-ang)   
    if 0<m:
        plt.plot([0,n*(ELT-ELF)],[ELF,ELF],'--',color='#000000',lw=2)
        plt.plot([0,-30],[ELF,ELF],'-',color='#000000',lw=1)
        plt.text(-30,ELF,'EL.{0:.3f}'.format(ELF),ha='left',va='bottom',fontsize=fsz)
        xs=-0.5*m*(ELF-ELB); ys=0.5*(ELF+ELB); ang=np.degrees(np.arctan(1/m))
        plt.text(xs-5,ys,'1:{0:.2f}'.format(m),va='center',ha='center',rotation=ang)   
    r=1.0
    theta=np.radians(np.arange(0,181,1))
    x=r*np.cos(theta)+xg-m*(ELF-ELB)
    y=r*np.sin(theta)+ELB+yg+r
    xx=np.append(x,[x[-1],x[0]])
    yy=np.append(y,[ELB+yg,ELB+yg])
    plt.fill(xx,yy,facecolor='#ffffff',edgecolor='#000000',lw=1)

    barrow(20,ELT,20,ELB)
    plt.text(20,0.5*(ELT+ELB),'H={0:.1f}m'.format(ELT-ELB),ha='right',va='center',fontsize=fsz,rotation=90)
    
    plt.plot([0,0],[ELB-5,ELT+15],'-.',color='#000000',lw=1)
    plt.text(0,ELT+15,'Dam axis',ha='center',va='bottom',fontsize=fsz,rotation=0)    
    
    plt.plot([tw,tw],[ELT,ELT+7],'-',color='#000000',lw=1)
    barrow(0,ELT+5,tw,ELT+5)
    plt.text(tw/2,ELT+7,'{0:.3f}'.format(tw),ha='center',va='bottom',fontsize=fsz,rotation=0)
    
    plt.plot([-m*(ELF-ELB),-m*(ELF-ELB)],[ELB,ELB-17],'-',color='#000000',lw=1)
    plt.plot([n*(ELT-ELB),n*(ELT-ELB)],[ELB,ELB-17],'-',color='#000000',lw=1)
    x1=-m*(ELF-ELB); x2=n*(ELT-ELB)
    barrow(x1,ELB-15,x2,ELB-15)
    plt.text(0.5*(x1+x2),ELB-15,'{0:.3f}'.format(x2-x1),ha='center',va='bottom',fontsize=fsz,rotation=0)

    xs=70; ys=ELT+15
    plt.text(xs,ys- 0,'UWL: U/S water level (FSL)',ha='left',va='center',fontsize=fsz)
    plt.text(xs,ys- 5,'USL: U/S sedimentation level',ha='left',va='center',fontsize=fsz)
    plt.text(xs,ys-10,'DWL: D/S water level',ha='left',va='center',fontsize=fsz)
    plt.text(xs,ys-15,'DSL: D/S sedimentation level',ha='left',va='center',fontsize=fsz)
    
    plt.tight_layout()
    plt.savefig(fnameF, dpi=100, bbox_inches="tight", pad_inches=0.1)
    plt.show()


def inp_BU():
    tw=7.5     # top width
    ELT=1850.0 # dam top level
    UWL=1845.0 # upstream water level
    USL=1805.8 # upstream sedimentation level
    DWL=1778.0 # downstream water level
    DSL=1778.0 # downstream sedimentation level
    ELB=1763.0 # dam foundation level
    ELF=ELB    # fillet level 
    n=0.85     # downstream gradient
    m=0.0     # fillet gradient
    xg=m*(ELF-ELB)+10.0 # location ogf gallery from upstream edge 
    yg=10.0    # height of gallery from dam foundation
    tstr='* Basochhu PSPP Upper dam non-overflow section (H={0:.1f}m)'.format(ELT-ELB)
    xmin=-60
    xmax=xmin+260
    dx=20
    ymin=1740
    ymax=ymin+140
    dy=10
    fpara=np.array([xmin,xmax,dx,ymin,ymax,dy])
    return tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m,fpara   


def inp_BL1():
    tw=7.5     # top width
    ELT=1120.0 # dam top level
    UWL=1115.0 # upstream water level
    USL=1078.1 # upstream sedimentation level
    DWL=1045.0 # dowstream water level
    DSL=1045.0 # downstream sedimentation level
    ELB=1025.0 # dam foundation level
    ELF=ELB   # fillet level 
    n=0.85     # downstream gradient
    m=0.00     # fillet gradient
    xg=m*(ELF-ELB)+10.0 # location ogf gallery from upstream edge 
    yg=10.0    # height of gallery from dam foundation
    if USL<ELB: USL=ELB
    if DWL<ELB: DWL=ELB
    if DSL<ELB: DSL=ELB
    tstr='* Basochhu PSPP Lower dam non-overflow section (H={0:.1f}m, Fillet elevation)'.format(ELT-ELB)
    xmin=-60
    xmax=xmin+260
    dx=20
    ymin=960
    ymax=ymin+190
    dy=10
    fpara=np.array([xmin,xmax,dx,ymin,ymax,dy])
    return tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m,fpara   


def inp_BL2():
    tw=7.5     # top width
    ELT=1120.0 # dam top level
    UWL=1115.0 # upstream water level
    USL=1078.1 # upstream sedimentation level
    DWL=1045.0 # dowstream water level
    DSL=1045.0 # downstream sedimentation level
    ELB=985.0 # dam foundation level
    ELF=1025.0    # fillet level 
    n=0.85     # downstream gradient
    m=0.20     # fillet gradient
    xg=m*(ELF-ELB)+10.0 # location ogf gallery from upstream edge 
    yg=10.0    # height of gallery from dam foundation
    tstr='* Basochhu PSPP Lower dam (non-overflow section (H={0:.1f}m)'.format(ELT-ELB)
    xmin=-60
    xmax=xmin+260
    dx=20
    ymin=960
    ymax=ymin+190
    dy=10
    fpara=np.array([xmin,xmax,dx,ymin,ymax,dy])
    return tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m,fpara   


def inp_JU():
    tw=7.5     # top width
    ELT=1345.0 # dam top level
    UWL=1340.0 # upstream water level
    USL=1305.4 # upstream sedimentation level
    DWL=1278.0 # downstream water level
    DSL=1278.0 # downstream sedimentation level
    ELB=1258.0 # dam foundation level
    ELF=ELB  # fillet level 
    n=0.85     # downstream gradient
    m=0.00     # fillet gradient
    xg=m*(ELF-ELB)+10.0    # location ogf gallery from upstream edge 
    yg=10.0    # height of gallery from dam foundation
    tstr='* Jerichhu PSPP Upper dam non-overflow section (H={0:.1f}m)'.format(ELT-ELB)
    xmin=-60
    xmax=xmin+260
    dx=20
    ymin=1240
    ymax=ymin+140
    dy=10
    fpara=np.array([xmin,xmax,dx,ymin,ymax,dy])
    return tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m,fpara   


def inp_JL1():
    tw=7.5     # top width
    ELT=675.0 # dam top level
    UWL=670.0 # upstream water level
    USL=620.4 # upstream sedimentation level
    DWL=540.0 # downstream water level
    DSL=540.0 # downstream sedimentation level
    ELB=560.0 # dam foundation level
    ELF=ELB  # fillet level 
    n=0.85     # downstream gradient
    m=0.00     # fillet gradient
    xg=m*(ELF-ELB)+10.0    # location ogf gallery from upstream edge 
    yg=10.0    # height of gallery from dam foundation
    if USL<ELB: USL=ELB
    if DWL<ELB: DWL=ELB
    if DSL<ELB: DSL=ELB
    tstr='* Jerichhu PSPP Lower dam non-overflow section (H={0:.1f}m, Fillet elevation)'.format(ELT-ELB)
    xmin=-60
    xmax=xmin+260
    dx=20
    ymin=500
    ymax=ymin+200
    dy=10
    fpara=np.array([xmin,xmax,dx,ymin,ymax,dy])
    return tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m,fpara   


def inp_JL2():
    tw=7.5     # top width
    ELT=675.0 # dam top level
    UWL=670.0 # upstream water level
    USL=620.4 # upstream sedimentation level
    DWL=540.0 # downstream water level
    DSL=540.0 # downstream sedimentation level
    ELB=520.0 # dam foundation level
    ELF=560.0  # fillet level 
    n=0.85     # downstream gradient
    m=0.50     # fillet gradient
    xg=m*(ELF-ELB)+10.0    # location ogf gallery from upstream edge 
    yg=10.0    # height of gallery from dam foundation
    tstr='* Jerichhu PSPP Lower dam non-overflow section (H={0:.1f}m)'.format(ELT-ELB)
    xmin=-60
    xmax=xmin+260
    dx=20
    ymin=500
    ymax=ymin+200
    dy=10
    fpara=np.array([xmin,xmax,dx,ymin,ymax,dy])
    return tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m,fpara   



def dambody(tw,ELT,ELB,ELF,n,m):
    # self weight of dam concrete
    d1=ELT-ELB
    d2=ELF-ELB
    w1=wc*tw*tw/n/2; x1=m*d2+tw/3*2; y1=d1-tw/n/3
    w2=wc*n*d1*d1/2; x2=m*d2+n*d1/3; y2=d1/3
    w3=wc*m*d2*d2/2; x3=m*d2/3*2   ; y3=d2/3
    vv=w1+w2+w3
    pp=vv*kh
    xx=(w1*x1+w2*x2+w3*x3)/vv
    yy=(w1*y1+w2*y2+w3*y3)/vv
    
    return vv,pp,xx,yy


def p_uwl(ELB,UWL):
    # U/S water pressure
    h=UWL-ELB
    vv=0.0
    pp=ww*h*h/2
    xx=0.0
    yy=h/3
    return vv,pp,xx,yy


def pd_uwl(ELB,UWL):
    # U/S dynamic water pressure
    h=UWL-ELB
    vv=0.0
    pp=ww*7/12*kh*h*h
    xx=0.0
    yy=0.4*h
    return vv,pp,xx,yy


def p_usl(ELB,USL):
    # U/S sedimentation pressure
    h=USL-ELB
    vv=0.0
    pp=ws*ce*h*h/2
    xx=0.0
    yy=h/3
    return vv,pp,xx,yy


def w_uwl(ELF,ELB,m,UWL):
    # U/S water weight
    if ELF<UWL:    
        h1=UWL-ELF
        h2=ELF-ELB
        w1=ww*h1*h2*m  ; x1=h2*m/2; y1=0.0
        w2=ww*h2*h2*m/2; x2=h2*m/3; y2=0.0
    else:
        h2=UWL-ELB
        w1=0.0         ; x1=0.0   ; y1=0.0
        w2=ww*h2*h2*m/2; x2=h2*m/3; y2=0.0
    vv=w1+w2
    pp=0.0
    xx=0.0
    yy=0.0
    if 0<m: 
        xx=(w1*x1+w2*x2)/vv
        yy=(w1*y1+w2*y2)/vv
    return vv,pp,xx,yy
    

def w_usl(ELF,ELB,m,USL):
    # U/S sedimentation weight
    if ELF<USL:    
        h1=USL-ELF
        h2=ELF-ELB
        w1=ws*h1*h2*m  ; x1=h2*m/2; y1=0.0
        w2=ws*h2*h2*m/2; x2=h2*m/3; y2=0.0
    else:
        h2=USL-ELB
        w1=0.0         ; x1=0.0   ; y1=0.0
        w2=ws*h2*h2*m/2; x2=h2*m/3; y2=0.0
    vv=w1+w2
    pp=0.0
    xx=0.0
    yy=0.0
    if 0<m:
        xx=(w1*x1+w2*x2)/vv
        yy=(w1*y1+w2*y2)/vv
    return vv,pp,xx,yy

    
def uplift(ELT,ELF,ELB,xg,yg,m,n,UWL,DWL):
    b=(ELF-ELB)*m+(ELT-ELB)*n
    h1=UWL-ELB
    h2=DWL-ELB
    hg=0.2*(h2-h1)+h2+yg
    u1=ww*xg*(h1-h2)/2    ; x1=xg/3       ; y1=0.0
    u2=ww*xg*h2           ; x2=xg/2       ; y2=0.0
    u3=ww*(b-xg)*(yg-h2)/2; x3=xg+(b-xg)/3; y3=0.0
    u4=ww*(b-xg)*h2       ; x4=xg+(b-xg)/2; y4=0.0
    uu=u1+u2+u3+u4
    pp=0.0
    xx=(u1*x1+u2*x2+u3*x3+u4*x4)/uu
    yy=(u1*y1+u2*y2+u3*y3+u4*y4)/uu
    return -uu,pp,xx,yy
    
    
def p_dwl(ELB,DWL):
    # D/S water pressure
    h=DWL-ELB
    vv=0.0
    pp=ww*h*h/2
    xx=0.0
    yy=h/3    
    return vv,-pp,xx,yy


def p_dsl(ELB,DSL):
    # D/S sedimentation pressure
    h=DSL-ELB
    vv=0.0
    pp=ws*ce*h*h/2
    xx=0.0
    yy=h/3    
    return vv,-pp,xx,yy


def w_dwl(ELT,ELF,ELB,n,m,DWL):
    # D/S water weight
    b=(ELF-ELB)*m+(ELT-ELB)*n
    h=DWL-ELB
    vv=ww*n*h*h/2
    pp=0.0
    xx=b-n*h/3
    yy=0.0
    return vv,pp,xx,yy


def w_dsl(ELT,ELF,ELB,n,m,DSL):
    # D/S sedimentation weight
    b=(ELF-ELB)*m+(ELT-ELB)*n
    h=DSL-ELB
    vv=ws*n*h*h/2
    pp=0.0
    xx=b-n*h/3
    yy=0.0
    return vv,pp,xx,yy


def calc(tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m):
    num=11
    fv=np.zeros(num,dtype=np.float64)
    fh=np.zeros(num,dtype=np.float64)
    x =np.zeros(num,dtype=np.float64)
    y =np.zeros(num,dtype=np.float64)
    si=['dam body weught',
        'U/S water pressure',
        'U/S dynamic pressure',
        'U/S sedimentation pressure',
        'U/S water weight',
        'U/S sedimentation weight',
        'Uplift',
        'D/S water pressure',
        'D/S sedimentation pressure',
        'D/S water weight',
        'D/S sedimentation weight'
    ]
    i=0;  fv[i],fh[i],x[i],y[i]=dambody(tw,ELT,ELB,ELF,n,m) # self-weight of concrete
    i=1;  fv[i],fh[i],x[i],y[i]=p_uwl(ELB,UWL)              # U/S water pressure
    i=2;  fv[i],fh[i],x[i],y[i]=pd_uwl(ELB,UWL)             # U/S dynamic water pressure
    i=3;  fv[i],fh[i],x[i],y[i]=p_usl(ELB,USL)              # U/S sedimentation pressure
    i=4;  fv[i],fh[i],x[i],y[i]=w_uwl(ELF,ELB,m,UWL)        # U/S water weight
    i=5;  fv[i],fh[i],x[i],y[i]=w_usl(ELF,ELB,m,USL)        # U/S sedimentation weight
    i=6;  fv[i],fh[i],x[i],y[i]=uplift(ELT,ELF,ELB,xg,yg,m,n,UWL,DWL) # uplift
    i=7;  fv[i],fh[i],x[i],y[i]=p_dwl(ELB,DWL)              # D/S water pressure
    i=8;  fv[i],fh[i],x[i],y[i]=p_dsl(ELB,DSL)              # D/S sedimentation pressure
    i=9;  fv[i],fh[i],x[i],y[i]=w_dwl(ELT,ELF,ELB,n,m,DWL)  # D/S water weight
    i=10; fv[i],fh[i],x[i],y[i]=w_dsl(ELT,ELF,ELB,n,m,DSL)  # D/S sedimentation weight
    fv=np.round(fv,decimals=3)
    fh=np.round(fh,decimals=3)
    x =np.round(x,decimals=3)
    y =np.round(y,decimals=3)
    fm=np.round(fv*x+fh*y,decimals=3)
    sum_v=np.sum(fv)
    sum_h=np.sum(fh)
    sum_m=np.sum(fm)
    bb=np.round(m*(ELF-ELB)+n*(ELT-ELB),decimals=3)
    ee=np.round(sum_m/sum_v-bb/2,decimals=3)
    pu=np.round(sum_v/bb*(1-6*ee/bb),decimals=3)
    pd=np.round(sum_v/bb*(1+6*ee/bb),decimals=3)
    sf=np.round((cc*bb+sum_v*ff)/sum_h,decimals=3)
    b6=np.round(bb/6,decimals=3)
    sj1='NG'
    sj2='NG'
    if 0<pu: sj1='ok'
    if 4<sf: sj2='ok'    

    ss=[]    
    ss=ss+[tstr]
    s='{0:30s},{1:>14s},{2:>14s},{3:>14s},{4:>14s},{5:>14s}'.format('Item','V(kN/m)','H(kN/m)','x(m)','y(m)','M(kN-m/m)')
    ss=ss+[s]
    for i in range(len(si)):
        s='{0:30s},{1:14.3f},{2:14.3f},{3:14.3f},{4:14.3f},{5:14.3f}'.format(si[i],fv[i],fh[i],x[i],y[i],fm[i])
        ss=ss+[s]
    s='{0:30s},{1:14.3f},{2:14.3f},{3:14s},{4:14s},{5:14.3f}'.format('Sum',sum_v,sum_h,'','',sum_m)
    ss=ss+[s]
    s='{0:30s},{1:>14s},{2:>14s},{3:>14s},{4:>14s},{5:>14s}'.format('Overturning','b/6(m)','e(m)','pu(kN/m2)','pd(kN/m2)','')
    ss=ss+[s]
    s='{0:30s},{1:14.3f},{2:14.3f},{3:14.3f},{4:14.3f},{5:>14s}'.format('',b6,ee,pu,pd,sj1)
    ss=ss+[s]
    s='{0:30s},{1:>14s},{2:>14s},{3:>14s},{4:>14s},{5:>14s}'.format('Sliding','c(kN/m2)','f','b(m)','SF','')
    ss=ss+[s]
    s='{0:30s},{1:14.3f},{2:14.3f},{3:14.3f},{4:14.3f},{5:>14s}'.format('',cc,ff,bb,sf,sj2)
    ss=ss+[s]
    ss=ss+['']
    for i in range(len(ss)):
        print(ss[i])

    return ss


def main():
    fnameF='fig_sec_BU.png'
    tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m,fpara=inp_BU()
    s_BU=calc(tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m)
    drawfig(tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m,fpara,fnameF)

    fnameF='fig_sec_BL.png'
    # above fillet
    tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m,fpara=inp_BL1()  
    s_BL1=calc(tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m)
    # below fillet
    tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m,fpara=inp_BL2()  
    s_BL2=calc(tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m)
    drawfig(tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m,fpara,fnameF)
    
    fnameF='fig_sec_JU.png'
    tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m,fpara=inp_JU()   
    s_JU=calc(tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m)
    drawfig(tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m,fpara,fnameF)

    fnameF='fig_sec_JL.png'
    # above fillet
    tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m,fpara=inp_JL1()   
    s_JL1=calc(tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m)
    #below fillet
    tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m,fpara=inp_JL2()   
    s_JL2=calc(tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m)
    drawfig(tstr,tw,ELT,UWL,USL,DWL,DSL,ELB,ELF,xg,yg,n,m,fpara,fnameF)

    fnameW='test.txt'
    f=open(fnameW,'w')
    for fn in [s_BU,s_BL1,s_BL2,s_JU,s_JL1,s_JL2]:
        for i in range(len(fn)):
            print(fn[i],file=f)
    f.close()    
    

#---------------
# Execute
#---------------
if __name__ == '__main__': main()

Sample of result is shown below.

* Basochhu PSPP Upper dam non-overflow section (H=87.0m)
Item                          ,       V(kN/m),       H(kN/m),          x(m),          y(m),     M(kN-m/m)
dam body weught               ,     76372.961,     11455.944,        24.450,        29.561,   2205968.057
U/S water pressure            ,         0.000,     33620.000,         0.000,        27.333,    918935.460
U/S dynamic pressure          ,         0.000,      5883.500,         0.000,        32.800,    192978.800
U/S sedimentation pressure    ,         0.000,      5495.520,         0.000,        14.267,     78404.584
U/S water weight              ,         0.000,         0.000,         0.000,         0.000,         0.000
U/S sedimentation weight      ,         0.000,         0.000,         0.000,         0.000,         0.000
Uplift                        ,    -12843.750,         0.000,        28.905,         0.000,   -371248.594
D/S water pressure            ,         0.000,     -1125.000,         0.000,         5.000,     -5625.000
D/S sedimentation pressure    ,         0.000,      -675.000,         0.000,         5.000,     -3375.000
D/S water weight              ,       956.250,         0.000,        69.700,         0.000,     66650.625
D/S sedimentation weight      ,      1147.500,         0.000,        69.700,         0.000,     79980.750
Sum                           ,     65632.961,     54654.964,              ,              ,   3162669.682
Overturning                   ,        b/6(m),          e(m),     pu(kN/m2),     pd(kN/m2),              
                              ,        12.325,        11.212,        80.148,      1694.915,            ok
Sliding                       ,      c(kN/m2),             f,          b(m),            SF,              
                              ,      3200.000,         1.200,        73.950,         5.771,            ok

Thank you.

Discharge capacity of spillway

Ishii-Fujimoto's formulas are shown below which is described in the Collection of Hydraulic Formulas (JSCE, 1985).

\begin{align} &Q=n\cdot C'\cdot B\cdot H^{3/2}\\ &C'=C_d\cdot\left\{1-M_d\cdot \left(\cfrac{H}{H_d}\right)^{1.5}\right\} \\ &C_d=1.971+0.498\cdot \xi+6.63\cdot \xi^2 \\ &M_d=0.0756\cdot\left(\cfrac{H_d}{B}\right)^{0.5}\cdot \left\{\cfrac{1}{n}+1.465\cdot \cfrac{n-1}{n}\cdot \left(\cfrac{b}{S'}\right)^{1.7}\right\} \end{align}

f:id:damyarou:20201025163415p:plain

\begin{align} & Q & & \text{dicharge} & & C_d & & \text{discharge coefficient at $H=H_d$} \\ & n & & \text{number of span} & & M_d & & \text{correction coefficient of effect of pier} \\ & C' & & \text{discharge coefficient} & & \xi & & \text{Value of $y/H_d$ at $x/H_d=0.5$} \\ & B & & \text{clear span} & & b & & \text{width of pier} \\ & H & & \text{overflow depth} & & S' & & \text{horizontal distance from top of overflow crest to upstream edge of pier} \\ & H_d & & \text{design head} & & & & \end{align}

For parabolic crest shape, the value of \xi can be obtained as follow.

\begin{equation} y=K\cdot x^2 \quad \Rightarrow \quad \cfrac{y}{H_d}=H_d\cdot K\cdot \left(\cfrac{x}{H_d}\right)^2 \quad \Rightarrow \quad \xi=\left(\cfrac{y}{H_d}\right)_{x/H_d=0.5}=H_d\cdot K\cdot \left(0.5\right)^2 \end{equation} 
import numpy as np


def calc_q():
    nsl=np.array([0.85, 0.85, 0.85, 0.85]) # downstream slope
    ELT=np.array([1850,1120,1345,675])     # dam top level
    FSL=np.array([1845,1115,1340,670])     # FSL
    hd=np.array([6.0, 19.0, 6.0, 19.0])    # design head
    ELC=FSL-hd                             # overflow crest level
    aa=ELT-ELC                             # dam top lebel - overflow crest level
    kk=1/4/aa/nsl**2                       # coefficient of ogee curve (parabola)
    nn=np.array([2, 6, 2, 10])             # number of gates
    bb=np.array([5.5, 12.0, 4.5, 12.0])    # clear span of pier
    bp=np.array([2.5, 4.0, 2.5, 4.0])      # pier width
    
    ss=bp/2+0.282*hd # distance from crest top to upstream edge of pier
    hh=hd
    md=0.0756*(hd/bb)**0.5*(1/nn+1.465*(nn-1)/nn*(bp/ss)**1.7)
    xi=hd*kk*0.5**2
    cd=1.971+0.498*xi+6.63*xi**2
    cc=cd*(1-md*(hh/hd)**1.5) # discharge coefficient
    qq=nn*cc*bb*hh**(3/2)     # discharge

    xxi=1/2/kk/nsl       # end point of ogee curve
    yyi=1/4/kk/nsl**2    # end point of ogee curve
    kki=1/kk             # inverse of palabola ogee curve coefficient 
    ctw=bp*(nn-1)+bb*nn   # cute width
    
    print('FSL ,{0:10.3f},{1:10.3f},{2:10.3f},{3:10.3f}'.format(FSL[0],FSL[1],FSL[2],FSL[3]))
    print('ELC ,{0:10.3f},{1:10.3f},{2:10.3f},{3:10.3f}'.format(ELC[0],ELC[1],ELC[2],ELC[3]))
    print('hd  ,{0:10.3f},{1:10.3f},{2:10.3f},{3:10.3f}'.format( hd[0], hd[1], hd[2], hd[3]))
    print('nn  ,{0:10.0f},{1:10.0f},{2:10.0f},{3:10.0f}'.format( nn[0], nn[1], nn[2], nn[3]))
    print('bb  ,{0:10.3f},{1:10.3f},{2:10.3f},{3:10.3f}'.format( bb[0], bb[1], bb[2], bb[3]))
    print('bp  ,{0:10.3f},{1:10.3f},{2:10.3f},{3:10.3f}'.format( bp[0], bp[1], bp[2], bp[3]))
    print('ss  ,{0:10.3f},{1:10.3f},{2:10.3f},{3:10.3f}'.format( ss[0], ss[1], ss[2], ss[3]))
    print('cc  ,{0:10.3f},{1:10.3f},{2:10.3f},{3:10.3f}'.format( cc[0], cc[1], cc[2], cc[3]))
    print('qq  ,{0:10.0f},{1:10.0f},{2:10.0f},{3:10.0f}'.format( qq[0], qq[1], qq[2], qq[3]))
    print('kki ,{0:10.3f},{1:10.3f},{2:10.3f},{3:10.3f}'.format(kki[0],kki[1],kki[2],kki[3]))
    print('ctw ,{0:10.3f},{1:10.3f},{2:10.3f},{3:10.3f}'.format(ctw[0],ctw[1],ctw[2],ctw[3]))
    
    return hd,ELC,xxi,yyi,kki

    
def main():
    hd,ELC,xxi,yyi,kki=calc_q()


#---------------
# Execute
#---------------
if __name__ == '__main__': main()

Results are shown below.

FSL ,  1845.000,  1115.000,  1340.000,   670.000
ELC ,  1839.000,  1096.000,  1334.000,   651.000
hd  ,     6.000,    19.000,     6.000,    19.000
nn  ,         2,         6,         2,        10
bb  ,     5.500,    12.000,     4.500,    12.000
bp  ,     2.500,     4.000,     2.500,     4.000
ss  ,     2.942,     7.358,     2.942,     7.358
cc  ,     1.842,     1.920,     1.824,     1.926
qq  ,       298,     11449,       241,     19143
kki ,    31.790,    69.360,    31.790,    69.360
ctw ,    13.500,    92.000,    11.500,   156.000

Thank you.

Hydraulic jump type energy dissipator

\begin{equation} v_1=0.9\cdot \sqrt{2\cdot g\cdot H} \end{equation}
\begin{equation} d_1=\cfrac{Q}{b\cdot v_1} \end{equation}
\begin{equation} F_1=\cfrac{v_1}{\sqrt{g\cdot d_1}} \end{equation}
\begin{equation} \cfrac{d_2}{d_1}=\cfrac{1}{2}\cdot \left(\sqrt{1+8\cdot F_1{}^2}-1\right) \end{equation}
\begin{equation} \cfrac{W}{d_1}=\cfrac{(1+2\cdot F_1{}^2)\cdot \sqrt{1+8\cdot F_1{}^2}-1-5\cdot F_1{}^2} {1+4\cdot F_1{}^2-\sqrt{1+8\cdot F_1{}^2}} -\left(\cfrac{\sqrt{g}}{C}\cdot F_1\right)^{2/3} \end{equation}
\begin{equation} L=4.5\cdot d_2 \end{equation}
\begin{align} & Q & & \text{discharge}\\ & b & & \text{width of apron}\\ & H & & \text{height from apron level to flood water level of reservoir}\\ & F_1 & & \text{Froude number at starting point of jump}\\ & v_1 & & \text{velocity at starting point of jump}\\ & d_1 & & \text{water depth at staring point of jump}\\ & d_2 & & \text{water depth at end point of jump}\\ & W & & \text{height of end sill}\\ & C & & \text{discharge coefficient of end sill}\\ & L & & \text{length of energy dissipator}\\ & g & & \text{gravity acceleration} \end{align} 
import numpy as np

def main():
    g=9.8
    q=42.0
    h0=554.1-321.0
    v1=26.5
    b=4.0

    d1=np.round(q/b/v1,decimals=3)
    fr=np.round(v1/np.sqrt(g*d1),decimals=3)
    d1d2=np.round(1/2*(np.sqrt(1+8*fr**2)-1),decimals=3)
    d2=np.round(d1*d1d2,decimals=3)
    print('v1=',v1)
    print('d1=',d1)
    print('Fr=',fr)
    print('d1d2=',d1d2)
    print('d2=',d2)

    c=1.704
    c1=(1+2*fr**2)*np.sqrt(1+8*fr**2)-1-5*fr**2
    c2=1+4*fr**2-np.sqrt(1+8*fr**2)
    cc=np.round(c1/c2-(np.sqrt(g)/c*fr)**(2/3),decimals=3)
    ww=np.round(d1*cc,decimals=3)
    print('cc=',cc)
    print('W=',ww)
    print('L=',4.5*d2)


#==============
# Execution
#==============
if __name__ == '__main__': main()

Calculation results are shown below.

v1= 26.5
d1= 0.396
Fr= 13.452
d1d2= 18.531
d2= 7.338
cc= 10.31
W= 4.083
L= 33.021

Thank you.

Normal depth of steep channel

\begin{gather} Q=A\cdot v \\ v=\cfrac{1}{n}\cdot R^{2/3}\cdot (\sin\theta)^{1/2} \\ A=b\cdot h \qquad R=\cfrac{b\cdot h}{b+2\cdot h} \end{gather}
Q(known) discharge
b(known) width of rectangular open channe
n(known) Manning's roughness coefficient
\theta(known) gradient of channel invert
h(unknown) normal depth

The critical depth of steep channel is calculated as following equation.

\begin{equation} h_c=\left(\cfrac{Q^2}{g\cdot b^2\cdot \cos\theta}\right)^{1/3} \end{equation} 
# normal depth and critical depth of rectangular cross section
import numpy as np
from scipy import optimize


def cal_hc(q,b,cs):
    # critical depth
    g=9.8
    hc=(q**2/g/b**2/cs)**(1/3)
    return hc


def func(h,q,b,n,sn):
    f=q-b*h/n*(b*h/(b+2*h))**(2/3)*sn**(1/2)    
    return f    


def main():
    q=42.0  # discharge
    b=4.0   # channel width
    n=0.014 # Manning's roughness coefficient
    theta=np.radians(37.0)
    sn=np.sin(theta)
    cs=np.cos(theta)
    
    h1=0.0
    h2=10.0
    hh=optimize.brentq(func,h1,h2,args=(q,b,n,sn))
    v=q/b/hh
    
    print('hn=',hh) # normal depth
    print('v=',v)
    
    hc=cal_hc(q,b,cs)
    print('hc=',hc) # critical depth
    
#==============
# Execution
#==============
if __name__ == '__main__': main()

Calculation results are shown below.

hn= 0.3962334595851377
v= 26.499528866122652
hc= 2.415097316241126

Thank you.