Programming

Programming

解説

値の変数への代入(:=)

Mapleは変数の初期設定で型宣言をする必要がない．数式処理の章で示したとおり，変数への代入は:=を使う．変数a,bにそれぞれ10,3を代入し，a+bの結果をc に代入するというプログラムは以下の通り．

> a:=10: b:=3: c:=a+b;

$$ c\, := \,13 $$

整数と浮動小数点数

浮動小数点数から整数に直すにはいくつかの関数がある.

• trunc：数値から数直線で 0 に向って最も近い整数
• round：数値の四捨五入
• floor：数値より小さな最も大きな整数
• ceil：数値より大きな最も小さな整数

負の値の時に floor と trunc は違った値を返す.

小数点以下を取りだすにはfrac が用意されている.

> frac(1.7);

$$ 0.7 $$

整数の割り算はirem(余り)とiquo(商).

> irem(7,3); #res: 1
> iquo(7,3); #res: 2

出力(print, printf)

Mapleではデフォルトで結果が出力される．これを抑えるには行末の”；”を”：”に変える必要がある．出力を明示的におこなうにはprintを使う．デバッグの時に便利．

> x:=1: print(x); #res: 1

さらに，出力を整えるのに便利なprintf関数がある．これはC言語と同じ構文で，

> printf("Hello world!!\n");

Hello world!!

と打ち込んでenterを押せば，出力が即座に表示される．値を表示するときには，

> i:=3: printf("%3d\n",i);

$$ 3 $$

となる．これは

「変数iに入っている値を，3桁の整数形式で打ち出した後，改行せよ」

と言う意味．%3dが出力の形式，\nが改行を意味する．OSによっては，\は¥と画面あるいはキーボードで表示されているかもしれない．実数の出力指定は%10.5fで，全部で10桁，小数点以下5桁で浮動小数点数を表示．複数の変数の出力は

> printf("%3d : %10.5f \n",i,a);

などとなる．

caption:printfの書式指定

%指定	意味
%o	整数を8進数で表示．
%d	整数を10進数で表示．
%x,%X	整数を16進数で表示．xは小文字，Xは大文字を使用．
%f	浮動小数点数として表示．
%e,%E	指数形式で表示．eは小文字，Eは大文字を使用．
%s	文字列を出力．

解説

for-loop

繰り返す操作はloopでおこなう．もっとも単純なfor-loop．

> for i from 1 to 3 do 
    i; 
  end do;

$$ 1 \notag \\ 2 \notag \\ 3 \notag $$

初期値や増減を調整したfor-loop

> for i from 10 by -2 to 0 do
    i; 
  end do;

$$ 10 \notag \\ 8 \notag \\ 6 \notag \\ 4 \notag \\ 2 \notag \\ 0 \notag $$

loop回数が少ないときは，loopの中身も出力される．これを止めるには，end do;の最後のセミコロンをコロンに変える．

二重ループ

i,jという二つの変数を使って2重化したループ．

> for i from 1 to 3 do 
    for j from 1 to 3 do 
      print(i,j); 
    end do; 
  end do;

$$ 1, 1 \notag \\ 1, 2 \notag \\ 1, 3 \notag \\ 2, 1 \notag \\ 2, 2 \notag \\ 2, 3 \notag \\ 3, 1 \notag \\ 3, 2 \notag \\ 3, 3 \notag $$

while-loop も同じように使える.

> i:=0; 
  while i<5 do 
    i:=i+1;
  end do;

$$ 0 \notag \\ 1 \notag \\ 2 \notag \\ 3 \notag \\ 4 \notag \\ 5 \notag $$

課題

1. printfを使って次のように表示せよ．

i) Hello world. ii) 1+1=2

1. 次の数を順に表示せよ．

i) 1から5までの整数．ii) 5から1までの整数 ．iii) 1から10にある偶数．

1. 9x9表を作れ．
2. 1 から 5 までの和を求めよ.
3. nを5にして，$n!=n \times (n-1) \times (n-1) \cdots 3 \times 2 \times 1$を求めよ.

解答例

> printf("Hello world!!\n");

Hello world!!

> i:=1; 
> printf("%d+%d=%d\n",i,i,i+i);

$$ 1 \notag \\ 1+1=2 \notag $$

i)

> for i from 1 to 5 do
    i;
  end do;

$$ 1 \notag \\ 2 \notag \\ 3 \notag \\ 4 \notag \\ 5 \notag $$

ii)

> for i from 5 to 1 by -1 do
    i; 
  end do;

$$ 5 \notag \\ 4 \notag \\ 3 \notag \\ 2 \notag \\ 1 \notag $$

iii)

> for i from 2 to 10 by 2 do
    i; 
  end do;

$$ 2 \notag \\ 4 \notag \\ 6 \notag \\ 8 \notag \\ 10 \notag $$

> for i from 1 to 9 do 
    for j from 1 to 9 do
      printf("%4d",i*j); 
    end do;
    printf("\n"); 
  end do;

   1   2   3   4   5   6   7   8   9
   2   4   6   8  10  12  14  16  18
   3   6   9  12  15  18  21  24  27
   4   8  12  16  20  24  28  32  36
   5  10  15  20  25  30  35  40  45
   6  12  18  24  30  36  42  48  54
   7  14  21  28  35  42  49  56  63
   8  16  24  32  40  48  56  64  72
   9  18  27  36  45  54  63  72  81

> sum1:=0; for i from 1 to 5 do 
    sum1:=sum1+i; 
  end do;

$$ 0 \notag \\ 1 \notag \\ 3 \notag \\ 6 \notag \\ 10 \notag \\ 15 \notag $$

> n:=5: 
  total1:=1:
  for i from 1 to n do 
    total1:=total1*i; 
  end do;

$$ 1 \notag \\ 2 \notag \\ 6 \notag \\ 24 \notag \\ 120 \notag $$

解説

配列は変数を入れる箱が沢山用意されていると考えればよい．配列を使うときは，箱を指す数(示数，index)をいじっているのか，箱の中身(要素)をいじっているのかを区別すれば，動作を理解しやすい．Mapleにはいくつかの配列構造が用意されている．もっとも，頻繁に使うlistを示す．

基本

リスト構造は，中に入れる要素を[]でくくる．

> restart; list1:=[1,3,5,7];

$$ {\it list1}\, := \,[1,3,5,7] $$

要素にアクセスするには，以下のようにインデックスを指定する．

> list1[2]; list1[-1]; list1[2..4];

$$ 3 \notag \\ 7 \notag \\ [3, 5, 7] \notag $$

-1,-2等は後ろから1番目，2番目を指す．C言語と違い0番目はない．

> list1[0];

Error, invalid subscript selector

ひとつの要素を書き換えるには，以下のようにする．

> list1[3]:=x: list1;

$$ [1, 3, x, 7] $$

要素の数，および要素の中身を取り出すには以下のようにする．

> nops(list1); 
> op(list1);

$$ 4 \notag \\ 1, 3, x, 7 \notag $$

for-loopの省略形

for-loopを省略するのによく使う手を二つ. (#より後ろはコメント文です)

配列の生成(seq)

> aa:=[]; #空で初期化 
  for i from 1 to 3 do 
    aa:=[op(aa),i]; #付け足していく
  end do:
  print(aa);

$$ {\it aa}\, := \,[] \notag \\ [1, 2, 3] \notag $$

同じことをseqを使って短く書くことができる．

> aa :=[seq(i,i=1..3)];

$$ {\it aa}\, := \,[1,2,3] $$

配列の和(sum)

> n:=nops(aa): 
  total:=0: 
  for i from 1 to n do 
    total:=total+aa[i]; 
  end do:
  print(total):

$$ 6 $$

同じことをsumを使って短く書くことができる．

> sum(aa[i],i=1..nops(aa));

Error, invalid subscript selector

sumやseqを使っていると，このようなエラーがよくでる．これは，for-loopをまわすときにiに値が代入されているため引っかかる．変数を換えるか，iを初期化すればよい．

> i;

$$ 4 $$

> sum(aa[j],j=1..nops(aa));

$$ 6 $$

リストへの付け足し(append, prepend)

opを用いると，リストに新たな要素を前後，あるいは途中に付け足すことができる．

> list1:=[op(list1),9];

$$ {\it list1}\, := \,[1,3,x,7,9] $$

2つの要素の入れ替え

要素の３，４番目の入れ替えは以下の通り．

> tmp:=list1[3]: 
  list1[3]:=list1[4]: 
  list1[4]:=tmp: 
  list1;

$$ [1, 3, 7, x, 9] $$

2次元配列(listlist)

[ ] を二重化することで 2 次元の配列を作ることも可能で,リストのリスト (listlist) と呼ばれる.

> l2:=[[1,2,3,4],[1,3,5,7]];

$$ {\it l2}\, := \,[[1,2,3,4],[1,3,5,7]] $$

要素へのアクセスは以下の通り．

> l2[2]; l2[2,3]; l2[2][3];

$$ [1, 3, 5, 7] \notag \\ 5 \notag \\ 5 \notag $$

listの表示(listplot)

listに入っている数値を視覚化するのにはlistplotが便利．

> la:=[1,2,3,4,3,2,1]; 
  with(plots): 
  listplot(la);

$$ [1, 2, 3, 4, 3, 2, 1] $$

課題

1. 1から100までの整数のうち5個をランダムに含んだ配列を生成せよ．

1から6までのランダムな数を生成する関数は，

> roll:=rand(1..6):

として作ることができる．実行は次の通り．

> seq(roll(),i=1..10);

$$ 5, 2, 5, 6, 2, 3, 4, 4, 6, 5 $$

1. さいころを１００回振って，出た目１から６が何回出たかを表示せよ．
2. コイン６枚を一度に投げて，表向きの枚数を数えるプログラムを書け．
3. 0から9までの整数5個から5桁の整数を作れ．（１桁目が０になっても気にするな）
4. 小数点以下8桁のそれぞれの桁の数を配列に格納せよ．8桁の少数は以下のようにして作られる．
5. 255以下の10進数をランダムに生成して，8桁の2進数へ変換せよ．

整数の割り算には irem(余り) と iquo(商) がある. 使用法は以下の通り．

> irem(7,3); #res: 1
> iquo(7,3); #res: 2

解答例

> roll:=rand(1..100): 
  [seq(roll(),i=1..5)];

$$ [27, 96, 17, 90, 34] $$

> roll:=rand(1..6):

$$ [0, 0, 0, 0, 0, 0] $$

> for i from 1 to 100 do 
    i1:=roll(); 
    A[i1]:=A[i1]+1; 
  end do:
  A;

$$ [16, 18, 21, 18, 18, 9] $$

> toss:=rand(0..1): 
  n:=6: 
  up:=0: 
  for i from 1 to n do 
    up:=up+toss(); 
  end do:
  up;

$$ 3 $$

> roll:=rand(0..9): 
  n:=5: 
  A:=[seq(roll(),i=1..n)];

$$ [5, 7, 3, 7, 6] $$

> sum1:=0:
  for i from 1 to n do 
    sum1:=sum1*10+A[i]; 
  end do: 
  sum1;

$$ 57376 $$

> restart; 
  n:=8:
  roll:=rand(10^(n-1)..10^n): 
  B:=evalf(roll()/10^n,8); 
  A:=[]:

$$ 0.19550684 $$

> B:=10*B; 
  for i from 1 to n do 
    A:=[op(A),floor(B)]; 
    B:=(B-A[i])*10; 
  end do: 
  A;

$$ 1.95506840 \notag \\ [1, 9, 5, 5, 0, 6, 8, 4] \notag $$

> n:=8: 
  roll:=rand(0..2^n-1): 
  B:=roll();

$$ 246 $$

> A:=[seq(0,j=1..n)]: 
  for i from 1 to n do 
    A[n-i+1]:=irem(B,2); 
    B:=iquo(B,2);
  end do: 
  A;

$$ [1, 1, 1, 1, 0, 1, 1, 0] $$

解説

if

もっとも簡単なif文の例．

> x:=-4:
  if (x<0) then
 　　y:=-x; 
  end if;

$$ 4 $$

例外付き．

> x:=3:
  if (x<0) then 
    y:=-x; 
  else 
    y:=x; 
  end if;

$$ 3 $$

2個の条件がある例

> x:=3:
  if (x<0) then 
    y:=-x; 
  elif (x>5) then 
    y:=x; 
  else
    y:=2*x; 
  end if;

$$ 6 $$

条件文に使える式と意味

関係演算子は<, <=, >, >=, =, <>で表記される．論理演算子にはand, or, xor, notがある．その他にもブール値を返す関数としてimplies, evalb, type などいくつかあり，条件分岐に使える．

caption:条件分岐のいくつかの例

xとyの値が一致	(x=y)
xとyの値が一致しない	(x<>y)
条件文を複数つなぐ	((x>0) and (x<4))
	((x<0) or (x>4))
	not (x=0)

nextとbreak

do-loopの途中で流れを変更するための命令．nextはdo-loop を一回スキップ．breakはそこで do-loop を一つ抜ける．以下のコードの出力結果を参照．

> for i from 1 to 5 do 
    if (i=3) then 
	  next; 
    end if; 
    print(i); 
  end do:

#res: 1 2 4 5

> for i from 1 to 5 do
    if (i=3) then
      break; 
    end if; 
    print(i); 
  end do:

#res: 1 2

課題

1. 西暦を代入したら，明治，大正，昭和，平成で答えてくれるプログラムを作成せよ．西暦1868, 1912,1926,1989年をそれぞれの元年とする．
2. 整数を代入したら，それ以下の素数をすべて表示するプログラムを作れ．素数かどうかの判定はMapleコマンドのisprimeを用いよ．
3. pが素数でp+2も素数のとき，これらは双子の素数と呼ばれる．10以上，100以下の双子の素数を全部見つけて出力せよ．
4. 素数判定を原理から実現せよ．

ある数nが素数かどうか(自分自身の数nと1以外の数で割りきれないかどうか)を判定せよ．割り算の余り(剰余)はiremで求まる．例えば

> residue:=irem(9,2);

として変数residue(余りの英語)をprintfしてみよ．番兵を置いておいて，n-1から2までの数でnを次々と割っていき，一度でも割り切れれば番兵にマークをつける．ループが終わった後に番兵のマークを見て素数(prime number)かどうかを判定すればよい．

1. うるう年かどうかを表示するプログラムをかけ．\\

うるう年は4で割り切れる数の年．ただし，100で割り切れる年はうるう年でなく，400で割り切れる年はうるう年．

1. ゴールドバッハの予想 \\

「6以上の偶数は二つの素数の和として表わされる」という予想を100以下の偶数について検証せよ．あらかじめ100までの素数をリストアップしておいてそのなかから組み合わせを探すと便利．

解答例

> year:=1890; 
  if year<1868 then printf("明治より前です．\n"); 
  elif year<1912 then printf("明治%d年です．\n",year-1868+1); 
  elif year<1926 then printf("大正%d年です．\n",year-1912+1); 
  elif year<1989 then printf("昭和%d年です．\n",year-1926+1); 
  elif year<2011 then printf("平成%d年です．\n",year-1989+1); 
  else printf("今年より後です．\n"); 
  end;

明治23年です．

> n:=10:
  for i from 1 to n do 
    if (isprime(i)) then 
      print(i); 
    end if; 
  end do;

#res: 2 3 5 7

> for i from 10 to 100-2 do 
    if (isprime(i) and isprime(i+2)) then 
      print(i,i+2);
    end if; 
  end do;

                                    11, 13
                                    17, 19
                                    29, 31
                                    41, 43
                                    59, 61
                                    71, 73

> n:=12: 
  banpei:=0:
  for i from 2 to n-1 do 
    residue:=irem(n,i); 
    # print(n,residue): 
    if residue=0 then 
      banpei:=1; 
      break; 
    end if; 
  end do: 
  if banpei=1 then 
    printf("%d is not prime number.\n",n); 
  else 
    printf("%d is prime number.\n",n); 
  end if;

12 is not prime number.

> year:=[2010,1984,2004,1800,1900,1600,2000]: 
  for i from 1 to nops(year) do 
    if (irem(year[i],400)=0) then 
      printf("%d is a leap year.\n",year[i]); 
    elif (irem(year[i],4)=0) and (irem(year[i],100)<>0) then
      printf("%d is a leap year.\n",year[i]); 
    else printf("%d is not a leap year.\n",year[i]); 
    end if;
  end do;

2010 is not a leap year. 
1984 is a leap year. 
2004 is a leap year. 
1800 is not a leap year. 
1900 is not a leap year. 
1600 is a leap year. 
2000 is a leap year.

別解

> for i from 1 to nops(year) do 
    if (irem(year[i],4)=0) and ((irem(year[i],100)<>0) or (irem(year[i],400)=0)) then 
      printf("%d is a leap year.\n",year[i]); 
    else 
      printf("%d is not a leap year.\n",year[i]); 
    end if;
  end do;

略

> prime1:=[]; 
  for i from 1 to 100 do 
    if isprime(i) then 
      prime1:=[op(prime1),i];
    end if; 
  end do; 
  prime1;

[2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 
73, 79, 83, 89, 97]

> nops(prime1);

$$ 25 $$

> for i from 6 to 100 by 2 do 
    for j1 from 1 to nops(prime1) do 
      for j2 from 1 to nops(prime1) do 
        if i=(prime1[j1]+prime1[j2]) then
          print(i,prime1[j1],prime1[j2]); 
          break;
        end if 
      end do; 
      if j2<=nops(prime1) then
        break; 
      end if; 
    end do; 
  end do;

6, 3, 3
8, 3, 5
10, 3, 7
中略
98, 19, 79
100, 3, 97

解説

基本

複雑な手続きや，何度も繰り返すルーチンはprocで作る． procは以下のようにして作る．

ユーザ関数名:=proc(仮引数)
  動作
end proc;

> test1:=proc(a) 
    print(a); 
  end proc:

procの呼び出しは，以下のようになる．

> test1(13);

$$ 13 $$

仮引数としてはどんな型(変数や配列)でもよい．複数指定するときにはコンマで区切る．仮引数をprocの中で変更することはできない．下で示すglobalで取り込むか，local変数にコピーして使う．

戻り値

procの戻り値はreturnで指定される．return文がないときは，最後の動作結果が戻り値となる．

> test2:=proc(a) 
    return a+1; 
  end proc:

> test2(13);

$$ 14 $$

グローバル(大域)，ローカル(局所)変数

procの内部だけで使われるのがlocal，外部を参照するのがglobal．global,localを省略してもMapleが適当に判断してくれるが，あまり信用せず，明示的に宣言した方が良い．宣言の仕方は以下の通り．

変数名:=proc(引数)
  local 変数,変数...;
  global 変数,変数...;
  動作の記述
end proc;

課題

1. 三角形の面積

底辺と高さを引数として，面積を返す関数areaを作れ．

1. MyIsprime

前章の課題４で求めた素数判定プログラムをprocにせよ．

1. ルートの距離

二つの位置座標 x1:=[0.0, 0.0] x2:=[1.0, 1.0] から距離を求めるMyDistance関数を作れ．

次に，4つの位置座標 x[1]=[0.0, 0.0] x[2]=[1.0, 1.0] x[3]=[1.0, 0.0] x[4]=[0.0, 1.0] を読み込んで，座標順に[1,2,3,4,1]と巡る距離を求めよ．

1. 最大数

ランダムな整数が格納されたリストを受け取り，そのリスト中の最大数を返す関数MyMaxを作れ．1から100までのランダムな整数のリストは次のようにして作れる．

> roll:=rand(1..100): 
  n:=50: 
  A:=[seq(roll(),i=1..n)];

$$ A\, := \,[45,96,6,98,59,44,100,38,69,27,96,17,90,34,18,52 \notag \\ ,56,43,83,25,90,93,60,93,14,50,47,8,46,44,9,77,59 \notag \\ ,16,1,70,77,39,92,71,67,78,51,53,12,19,63,40,90,3 \notag $$

解答例

> area:=proc(base,height) 
    base*height/2; 
  end proc:

> area(3,4); #res: 6

> restart; 
  n:=19: 
  banpei:=0; 
  for i from 2 to n-1 do 
    amari:=irem(n,i);
    print(amari): 
    if amari=0 then 
      banpei:=1; 
      break;
    end if; 
  end do: 
  if banpei=1 then 
    print(n," is not prime number."); 
  else 
    print(n," is prime number."); 
  end if;

$$ 0 \notag \\ 1 \notag \\ 1 \notag \\ 1 \notag \\ 19, "\, is\, prime\, number." \notag $$

> MyIsprime:=proc(n) 
    local i,amari; 
    for i from 2 to evalf(sqrt(n)) do
      amari:=irem(n,i); 
      if amari=0 then 
        return false; 
      end if; 
    end do:
    return true;
  end proc:

> MyIsprime(104729);

$$ true $$

> restart; x1:=[0.0, 0.0]: x2:=[1.0, 1.0]:

> MyDistance:=proc(x1,x2) 
    local dx,dy; 
    dx:=(x1[1]-x2[1]); 
    dy:=(x1[2]-x2[2]);
    sqrt(dx^2+dy^2); 
  end proc:

> MyDistance(x1,x2);

$$ 1.414213562 $$

> x[1]:=[0.0, 0.0]: x[2]:=[1.0, 1.0]: x[3]:=[1.0, 0.0]: x[4]:=[0.0, 1.0]: x[5]:=x[1]: 
  sum(MyDistance(x[i],x[i+1]),i=1..4);

$$ 4.828427124 $$

> MyMax:=proc(A) 
    local imax,i; 
    imax:=0; 
    for i from 1 to nops(A) do 
      if A[i]>imax then 
          imax:=A[i];
      end if 
    end do; 
    return imax; 
  end proc:

> MyMax(A);

$$ 100 $$