teacup. [ 掲示板 ] [ 掲示板作成 ] [ 有料掲示板 ] [ ブログ ]

【重要】teacup. byGMOのサービス終了について
 <思いやりのあるコミュニティ宣言>
 teacup.掲示板は、皆様の権利を守りながら、思いやり、温かみのあるコミュニティづくりを応援します。
 いつもご協力いただきありがとうございます。

 投稿者
  題名
  内容 入力補助画像・ファイル<IMG>タグが利用可能です。(詳細)
    
 URL
[ ケータイで使う ] [ BBSティッカー ] [ 書込み通知 ] [ 検索 ]

スレッド一覧

  1. スレッドが使えます(2)
  2. Amusement_Program(12)
  3. Chromebook(3)
  4. 人の色覚の数理(14)
  5. Apple Silicon (M1プロセッサ)での動作(7)
  6. Full BASIC互換ライブラリ(7)
  7. Paract BASIC(22)
  8. 改修予定のJIS非互換(3)
  9. 複数ページ長編プログラム(新規投稿)(16)
  10. 「十進BASIC第2掲示板」投稿記事リスト(17)
  11. 「十進BASIC掲示板過去ログ」インデックス(トピック)(17)
  12. 「十進BASIC掲示板過去ログ」インデックス(ツリー)(91)
スレッド一覧(全12)  他のスレッドを探す  スレッド作成

*掲示板をお持ちでない方へ、まずは掲示板を作成しましょう。無料掲示板作成


第2掲示板過去ログ公開

 投稿者:SHIRAISHI Kazuo  投稿日:2022年 3月10日(木)11時50分39秒
返信・引用
  第2掲示板の過去ログを公開しました。
https://decimalbasic.ninja-web.net/bbs2/logs2.html
画像によってはクリックすると浮き出てくるものがあるのですが,それらはクリックしたときteacupのサーバからダウンロードされるようです。念のため,ダウンロードして保存してあります。万一のときにはそれをもとにリンクを貼り直します。
 
 

Re: 十進BASIC第3掲示板テスト

 投稿者:SHIRAISHI Kazuo  投稿日:2022年 3月 7日(月)10時04分26秒
返信・引用 編集済
  > No.5029[元記事へ]

> FC2に新掲示板を作成しました。
> 投稿テストへの協力をお願いします。
> https://deimal-basic.bbs.fc2.com/
> 英数字だけの本文は不可のようです。

FC2の掲示板に移行します。
第2掲示板はteacupのサービス終了まで存続させます。
過去ログを
https://decimalbasic.ninja-web.net/BBS2/logs2.html
に置きます。
 

Re: 十進BASIC第3掲示板テスト

 投稿者:SHIRAISHI Kazuo  投稿日:2022年 3月 7日(月)10時01分51秒
返信・引用 編集済
  > No.5028[元記事へ]

ZawaZawaの掲示板は凍結しました。
インデントが再現されるのはteacup同様で好ましいのですが,
不等号がタグと誤認されたり,*と*で囲むと*が消えたりなど扱いにくい点が見つかりました。
書き込みテストへのご協力ありがとうございました。

> 十進BASIC第3掲示板のテスト運用を開始しました。
> 投稿テストへの協力をお願いします。
> 150行制限があります。
>
> https://zawazawa.jp/decimal_basic/
 

十進BASIC第3掲示板テスト

 投稿者:SHIRAISHI Kazuo  投稿日:2022年 3月 6日(日)10時16分25秒
返信・引用 編集済
  FC2に新掲示板を作成しました。
投稿テストへの協力をお願いします。
https://deimal-basic.bbs.fc2.com/
英数字だけの本文は不可のようです。
 

十進BASIC第3掲示板テスト

 投稿者:SHIRAISHI Kazuo  投稿日:2022年 3月 5日(土)07時47分14秒
返信・引用 編集済
  十進BASIC第3掲示板のテスト運用を開始しました。
投稿テストへの協力をお願いします。
150行制限があります。

https://zawazawa.jp/decimal_basic/
 

Re: インストーラ版に同梱されていないファイルがあります

 投稿者:SHIRAISHI Kazuo  投稿日:2022年 2月28日(月)16時34分18秒
返信・引用
  > No.5026[元記事へ]

ご報告ありがとうございました。
修正版を作成します。


> Ver.7.8.6.4 の Windowsインストーラ版に下記のファイルが同梱されていません。
>
> テキストファイル  REVISION.TXT
> ファイルフォルダー  Tutorial
>
> アーカイブ版にはあります。
 

インストーラ版に同梱されていないファイルがあります

 投稿者:nagram  投稿日:2022年 2月28日(月)11時19分1秒
返信・引用
  Ver.7.8.6.4 の Windowsインストーラ版に下記のファイルが同梱されていません。

テキストファイル  REVISION.TXT
ファイルフォルダー  Tutorial

アーカイブ版にはあります。
 

フラグメントシェーダー 青海波模様

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時36分3秒
返信・引用
  フラグメントシェーダー  青海波模様

青海波模様です。


DIM ST(2),C(3),T1(2),T2(2),T3(2)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
CALL VEC2(T1,1,-.5)
CALL VEC2(T2,0,0)
CALL VEC2(T3,1,.5)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(ST,X,Y)
      MAT ST=5*ST
      FOR I=1 TO 2
         LET ST(I)=FRACT(ST(I))
      NEXT I
      MAT ST=2*ST
      LET INDEX=0
      LET INDEX=INDEX+STEP(1,MOD(ST(1),2))
      LET INDEX=INDEX+STEP(1,MOD(ST(2),2))*2
      FOR I=1 TO 2
         LET ST(I)=FRACT(ST(I))
      NEXT I
      IF INDEX=1 OR INDEX=3 THEN
         LET ST(1)=1-ST(1)
      END IF
      LET D=DISTANCE(T1,ST)
      LET T=0
      IF D<1 THEN
         LET T=ABS(COS(D*2.65*PI))
      END IF
      LET D=DISTANCE(T2,ST)
      IF D<1 AND T=0 THEN
         LET T=ABS(COS(D*2.65*PI))
      END IF
      LET D=DISTANCE(T3,ST)
      IF T=0 THEN
         LET T=ABS(COS(D*2.65*PI))
      END IF
      IF T>.5 THEN
         CALL VEC3(C,0.,.5,1)
      ELSE
         CALL VEC3(C,1,1,1)
      END IF
      CALL SETCOLOR(C)
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  SUB SETCOLOR(COL())
SET COLOR COLORINDEX(CLAMP(COL(1),0,1),CLAMP(COL(2),0,1),CLAMP(COL(3),0,1))
END SUB

EXTERNAL  FUNCTION DISTANCE(A(),B())
LET DISTANCE=SQR((B(1)-A(1))^2+(B(2)-A(2))^2)
END FUNCTION

EXTERNAL  FUNCTION FRACT(X)
LET FRACT=X-INT(X)
END FUNCTION

EXTERNAL  FUNCTION STEP(A,X)
IF X0 THEN
   LET SIGN=1
ELSE
   LET SIGN=0
END IF
END FUNCTION
 

フラグメントシェーダー 七宝模様

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時35分22秒
返信・引用
  フラグメントシェーダー  七宝模様

七宝模様です。


DIM ST(2),C(3),T1(2),T2(2)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
CALL VEC2(T1,0,0)
CALL VEC2(T2,1,1)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(ST,X,Y)
      MAT ST=5*ST
      FOR I=1 TO 2
         LET ST(I)=FRACT(ST(I))
      NEXT I
      MAT ST=2*ST
      LET INDEX=0
      LET INDEX=INDEX+STEP(1,MOD(ST(1),2))
      LET INDEX=INDEX+STEP(1,MOD(ST(2),2))*2
      FOR I=1 TO 2
         LET ST(I)=FRACT(ST(I))
      NEXT I
      IF INDEX=1 OR INDEX=2 THEN
         LET ST(1)=1-ST(1)
      END IF
      IF DISTANCE(ST,T1)<1 AND DISTANCE(ST,T2)<1 THEN
         CALL VEC3(C,1,1,0)
      ELSE
         CALL VEC3(C,1,0,0)
      END IF
      CALL SETCOLOR(C(1),C(2),C(3))
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  SUB SETCOLOR(R,G,B)
SET COLOR COLORINDEX(CLAMP(R,0,1),CLAMP(G,0,1),CLAMP(B,0,1))
END SUB

EXTERNAL  FUNCTION DISTANCE(A(),B())
LET DISTANCE=SQR((B(1)-A(1))^2+(B(2)-A(2))^2)
END FUNCTION

EXTERNAL  FUNCTION FRACT(X)
LET FRACT=X-INT(X)
END FUNCTION

EXTERNAL  FUNCTION STEP(A,X)
IF X0 THEN
   LET SIGN=1
ELSE
   LET SIGN=0
END IF
END FUNCTION
 

フラグメントシェーダー 市松模様

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時34分43秒
返信・引用
  フラグメントシェーダー  市松模様

市松模様を描いています。
https://junk-box.net/kuyo/index.php/2020/11/29/material_part1/

※「鬼滅の刃」竈門炭治郎のイメージ ?

DIM ST(2),C(3)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(ST,X,Y)
      !! CALL ROTATE(ST,RAD(25)) !注釈を外すと傾けることができます。
      MAT ST=5*ST
      FOR I=1 TO 2
         LET ST(I)=FRACT(ST(I))
      NEXT I
      MAT ST=2*ST
      LET INDEX=0
      LET INDEX=INDEX+STEP(1,MOD(ST(1),2))
      LET INDEX=INDEX+STEP(1,MOD(ST(2),2))*2
      IF INDEX=1 OR INDEX=2 THEN
         CALL VEC3(C,0,1,0)
      ELSE
         CALL VEC3(C,0,0,0)
      END IF
      CALL SETCOLOR(C(1),C(2),C(3))
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  SUB ROTATE(A(),T)
DIM M(2,2)
LET S=SIN(T)
LET C=COS(T)
LET M(1,1)=C
LET M(1,2)=S
LET M(2,1)=-S
LET M(2,2)=C
MAT A=A*M
END SUB

EXTERNAL  SUB SETCOLOR(R,G,B)
SET COLOR COLORINDEX(CLAMP(R,0,1),CLAMP(G,0,1),CLAMP(B,0,1))
END SUB

EXTERNAL  FUNCTION FRACT(X)
LET FRACT=X-INT(X)
END FUNCTION

EXTERNAL  FUNCTION STEP(A,X)
IF X0 THEN
   LET SIGN=1
ELSE
   LET SIGN=0
END IF
END FUNCTION
-----------------------------------------------------------------------------------------------------------
円形にしてみた


DIM P(2),COL(3)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
LET N=20
LET SIZE=2/N
LET NN=2
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(P,X,Y)
      LET D=NLENGTH(P,NN)
      IF X=0 AND Y=0 THEN LET TH=0 ELSE LET TH=ANGLE(X,Y)
      LET D=INT(D/SIZE)
      LET TH=INT(TH/SIZE/PI)
      IF MOD(D+TH,2)=0 THEN
         CALL VEC3(COL,0,0,0)
      ELSE
         CALL VEC3(COL,0,1,0)
      END IF
      CALL SETCOLOR(COL(1),COL(2),COL(3))
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  SUB SETCOLOR(R,G,B)
SET COLOR COLORINDEX(CLAMP(R,0,1),CLAMP(G,0,1),CLAMP(B,0,1))
END SUB

EXTERNAL  SUB VEC2(A(),X,Y)
LET A(1)=X
LET A(2)=Y
END SUB

EXTERNAL  SUB VEC3(A(),X,Y,Z)
LET A(1)=X
LET A(2)=Y
LET A(3)=Z
END SUB

EXTERNAL  FUNCTION CLAMP(X,A,B)
LET CLAMP=MIN(B,MAX(X,A))
END FUNCTION

EXTERNAL  FUNCTION NLENGTH(P(),N)
LET NLENGTH=(ABS(P(1))^N+ABS(P(2))^N)^(1/N)
END FUNCTION
 

フラグメントシェーダー パーリンノイズ

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時33分34秒
返信・引用
  フラグメントシェーダー  パーリンノイズ

ノイズ関数を定義してパーリンノイズを生成しています。
雲(霧?)のようなものを描きます。


DIM P(2),COL(3)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(P,X,Y)
      MAT P=3*P
      LET N=FBM(P)
      CALL VEC3(COL,N,N,N)
      CALL SETCOLOR(COL(1),COL(2),COL(3))
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  SUB SETCOLOR(R,G,B)
SET COLOR COLORINDEX(CLAMP(R,0,1),CLAMP(G,0,1),CLAMP(B,0,1))
END SUB

EXTERNAL  FUNCTION RANDOM(S())
DIM V(2)
CALL VEC2(V,12.9898,78.233)
LET RANDOM=FRACT(SIN(DOT(S,V))*43758.5453123)
END FUNCTION

EXTERNAL  FUNCTION NOISE(ST())
DIM I(2),F(2),U(2),I1(2),I2(2),I3(2)
FOR J=1 TO 2
   LET I(J)=FLOOR(ST(J))
   LET F(J)=FRACT(ST(J))
NEXT J
CALL VEC2(I1,1,0)
CALL VEC2(I2,0,1)
CALL VEC2(I3,1,1)
MAT I1=I1+I
MAT I2=I2+I
MAT I3=I3+I
LET A=RANDOM(I)
LET B=RANDOM(I1)
LET C=RANDOM(I2)
LET D=RANDOM(I3)
FOR J=1 TO 2
   LET U(J)=F(J)*F(J)*(3-2*F(J))
NEXT J
LET NOISE=MIX(A,B,U(1))+(C-A)*U(2)*(1-U(1))+(D-B)*U(1)*U(2)
END FUNCTION

EXTERNAL  FUNCTION FBM(ST())
LET VALUE=0
LET AMPLITUDE=.5
LET FREQUENCY=0
FOR I=0 TO 5
   LET VALUE=VALUE+AMPLITUDE*NOISE(ST)
   MAT ST=2*ST
   LET AMPLITUDE=AMPLITUDE*.5
NEXT I
LET FBM=VALUE
END FUNCTION

EXTERNAL  FUNCTION FRACT(X)
LET FRACT=X-INT(X)
END FUNCTION

EXTERNAL  FUNCTION FLOOR(X)
LET FLOOR=INT(X)
END FUNCTION

EXTERNAL  FUNCTION STEP(A,X)
IF X<A THEN LET STEP=0 ELSE LET STEP=1
END FUNCTION

EXTERNAL  SUB VEC2(A(),X,Y)
LET A(1)=X
LET A(2)=Y
END SUB

EXTERNAL  SUB VEC3(A(),X,Y,Z)
LET A(1)=X
LET A(2)=Y
LET A(3)=Z
END SUB

EXTERNAL  FUNCTION CLAMP(X,A,B)
LET CLAMP=MIN(B,MAX(X,A))
END FUNCTION

EXTERNAL  FUNCTION MIX(X,Y,A)
LET MIX=X*(1-A)+Y*A
END FUNCTION
 

フラグメントシェーダー ホワイトノイズ

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時32分48秒
返信・引用
  フラグメントシェーダー  ホワイトノイズ

乱数randomを定義してホワイトノイズ(砂嵐)を描いています。

randomはデダラメな数値を返す関数です。
引数が同じなら同じ値を返します。

https://nogson2.hatenablog.com/entry/2017/11/18/150645


DIM P(2),COL(3)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(P,X,Y)
      LET N=RANDOM(P)
      CALL VEC3(COL,N,N,N)
      CALL SETCOLOR(COL(1),COL(2),COL(3))
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  SUB SETCOLOR(R,G,B)
SET COLOR COLORINDEX(CLAMP(R,0,1),CLAMP(G,0,1),CLAMP(B,0,1))
END SUB

EXTERNAL  FUNCTION RANDOM(S()) !乱数
DIM V(2)
CALL VEC2(V,12.9898,78.233)
LET RANDOM=FRACT(SIN(DOT(S,V))*43758.5453123)
END FUNCTION

EXTERNAL  FUNCTION FRACT(X)
LET FRACT=X-INT(X)
END FUNCTION

EXTERNAL  SUB VEC2(A(),X,Y)
LET A(1)=X
LET A(2)=Y
END SUB

EXTERNAL  SUB VEC3(A(),X,Y,Z)
LET A(1)=X
LET A(2)=Y
LET A(3)=Z
END SUB

EXTERNAL  FUNCTION CLAMP(X,A,B)
LET CLAMP=MIN(B,MAX(X,A))
END FUNCTION
 

フラグメントシェーダー 三項ブーリアン演算

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時31分47秒
返信・引用
  フラグメントシェーダー  三項ブーリアン演算

ブーリアン演算を3項でやってみた。


DIM ST(2),C(3)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(ST,X,Y)
      IF DISTANCE(ST)<0 THEN
         CALL VEC3(C,1,0,0)
      ELSE
         CALL VEC3(C,0,0,0)
      END IF
      CALL SETCOLOR(C)
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  FUNCTION MIN3(X,Y,Z)
LET MIN3=MIN(X,MIN(Y,Z))
END FUNCTION

EXTERNAL  FUNCTION MIN4(X,Y,Z,W)
LET MIN4=MIN(MIN(X,Y),MIN(Z,W))
END FUNCTION

EXTERNAL  FUNCTION MAX3(X,Y,Z)
LET MAX3=MAX(X,MAX(Y,Z))
END FUNCTION

EXTERNAL  FUNCTION CIRCLE(P(),B(),R)
DIM PP(2)
MAT PP=P-B
LET CIRCLE=LENGTH(PP)-R
END FUNCTION

EXTERNAL  FUNCTION DISTANCE(P())
DIM P1(2),P2(2),P3(2)
CALL VEC2(P1,.3*COS(0),.3*SIN(0))
CALL VEC2(P2,.3*COS(2*PI/3),.3*SIN(2*PI/3))
CALL VEC2(P3,.3*COS(4*PI/3),.3*SIN(4*PI/3))
LET A=CIRCLE(P,P1,.5)
LET B=CIRCLE(P,P2,.5)
LET C=CIRCLE(P,P3,.5)
LET DISTANCE=MIN3(A,B,C)
!LET DISTANCE=MAX3(A,B,C)
!LET DISTANCE=MAX3(-A,B,C)
!LET DISTANCE=MAX3(A,-B,C)
!LET DISTANCE=MAX3(A,B,-C)
!LET DISTANCE=MAX3(A,-B,-C)
!LET DISTANCE=MAX3(-A,B,-C)
!LET DISTANCE=MAX3(-A,-B,C)
!LET DISTANCE=MIN3(MIN(MAX(-A,B),MAX(A,-B)),MIN(MAX(-B,C),MAX(B,-C)),MIN(MAX(-A,C),MAX(A,-C)))
!LET DISTANCE=MAX3(MIN(A,B),MIN(B,C),MIN(A,C))
!LET DISTANCE=MIN3(MAX(A,B),MAX(B,C),MAX(A,C))
!LET DISTANCE=MIN3(MAX(-A,B),MAX(-B,C),MAX(A,-C))
!LET DISTANCE=MIN3(MAX3(-A,-B,C),MAX3(A,-B,-C),MAX3(-A,B,-C))
!LET DISTANCE=MIN3(MAX3(-A,B,C),MAX3(A,-B,C),MAX3(A,B,-C))
!LET DISTANCE=MIN(MAX3(-A,B,C),MAX3(A,-B,C))
!LET DISTANCE=MIN3(MAX3(-A,-B,C),MAX3(A,B,C),MAX3(A,B,-C))
!LET DISTANCE=MIN3(MAX3(-A,B,C),MAX3(-A,-B,C),MAX3(A,-B,-C))
!LET DISTANCE=MIN3(MAX3(A,B,C),MAX3(-A,B,C),MAX3(A,-B,C))
!LET DISTANCE=MAX(A,MIN(-B,C))
!LET DISTANCE=MIN3(MAX3(A,-B,-C),MAX3(-A,B,-C),MAX3(A,B,C))
!LET DISTANCE=MIN3(MAX3(A,-B,-C),MAX3(-A,B,C),MAX3(A,B,C))
!LET DISTANCE=MIN4(MAX3(A,-B,C),MAX3(-A,-B,C),MAX3(A,B,-C),MAX3(A,-B,-C))
!LET DISTANCE=MIN3(MAX3(A,B,-C),MAX3(A,-B,-C),MAX3(-A,-B,C))
END FUNCTION

EXTERNAL  SUB SETCOLOR(COL())
SET COLOR COLORINDEX(CLAMP(COL(1),0,1),CLAMP(COL(2),0,1),CLAMP(COL(3),0,1))
END SUB

EXTERNAL  SUB VEC2(A(),X,Y)
LET A(1)=X
LET A(2)=Y
END SUB

EXTERNAL  SUB VEC3(A(),X,Y,Z)
LET A(1)=X
LET A(2)=Y
LET A(3)=Z
END SUB

EXTERNAL  FUNCTION CLAMP(X,A,B)
LET CLAMP=MIN(B,MAX(X,A))
END FUNCTION

EXTERNAL  FUNCTION LENGTH(A())
LET LENGTH=SQR(A(1)^2+A(2)^2)
END FUNCTION
 

フラグメントシェーダー ブーリアン演算 smooth min関数

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時30分46秒
返信・引用
  フラグメントシェーダー  ブーリアン演算 smooth min関数

ブーリアン演算のmin関数をsmooth min関数にすると
接合部をなめらかにすることができます。min関数と比べてみてください。
https://www.iquilezles.org/www/articles/smin/smin.htm

DIM ST(2),C(3)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(ST,X,Y)
      !IF ABS(DISTANCE(ST))<.01 THEN
      IF DISTANCE(ST)<0 THEN
         CALL VEC3(C,1,0,0)
      ELSE
         CALL VEC3(C,0,0,0)
      END IF
      CALL SETCOLOR(C)
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  FUNCTION CIRCLE(P(),B(),R)
DIM PP(2)
MAT PP=P-B
LET CIRCLE=LENGTH(PP)-R
END FUNCTION

EXTERNAL  FUNCTION DISTANCE(P())
DIM P1(2),P2(2)
CALL VEC2(P1,-.25,0)
CALL VEC2(P2,.25,0)
LET DISTANCE=SMIN(CIRCLE(P,P1,.6),CIRCLE(P,P2,.4))
!LET DISTANCE=SMIN2(CIRCLE(P,P1,.6),CIRCLE(P,P2,.4))
!LET DISTANCE=SMIN3(CIRCLE(P,P1,.6),CIRCLE(P,P2,.4))
!LET DISTANCE=SMIN4(CIRCLE(P,P1,.6),CIRCLE(P,P2,.4))
!LET DISTANCE=SMIN5(CIRCLE(P,P1,.6),CIRCLE(P,P2,.4))
END FUNCTION

EXTERNAL  FUNCTION SMIN(A,B) !exponential smooth min
LET K=32
LET RES=EXP(-K*A)+EXP(-K*B)
LET SMIN=-LOG(RES)/K
END FUNCTION

EXTERNAL  FUNCTION SMIN2(A,B) !polynomial smooth min
LET K=.1
LET H=CLAMP(.5+.5*(B-A)/K,0,1)
LET SMIN2=MIX(B,A,H)-K*H*(1-H)
END FUNCTION

EXTERNAL  FUNCTION SMIN3(A,B) !polynomial smooth min
LET K=.1
LET H=MAX(K-ABS(A-B),0)/K
LET SMIN3=MIN(A,B)-H*H*K/4
END FUNCTION

EXTERNAL  FUNCTION SMIN4(A,B) !polynomial smooth min
LET K=.1
LET H=MAX(K-ABS(A-B),0)/K
LET SMIN4=MIN(A,B)-H*H*H*K/6
END FUNCTION

EXTERNAL  FUNCTION SMIN5(A,B) !root smooth min
LET K=.01
LET H=A-B
LET SMIN5=.5*((A+B)-SQR(H*H+K))
END FUNCTION

!EXTERNAL  FUNCTION SMAX(A,B)
!LET SMAX=A+B-SMIN(A,B)
!END FUNCTION

EXTERNAL  SUB SETCOLOR(COL())
SET COLOR COLORINDEX(CLAMP(COL(1),0,1),CLAMP(COL(2),0,1),CLAMP(COL(3),0,1))
END SUB

EXTERNAL  SUB VEC2(A(),X,Y)
LET A(1)=X
LET A(2)=Y
END SUB

EXTERNAL  SUB VEC3(A(),X,Y,Z)
LET A(1)=X
LET A(2)=Y
LET A(3)=Z
END SUB

EXTERNAL  FUNCTION CLAMP(X,A,B)
LET CLAMP=MIN(B,MAX(X,A))
END FUNCTION

EXTERNAL  FUNCTION LENGTH(A())
LET LENGTH=SQR(A(1)^2+A(2)^2)
END FUNCTION

EXTERNAL  FUNCTION MIX(X,Y,A)
LET MIX=X*(1-A)+Y*A
END FUNCTION
 

フラグメントシェーダー ブーリアン演算

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時29分43秒
返信・引用
  フラグメントシェーダー ブーリアン演算

円を2つ定義してブーリアン演算により合体させています。
ブーリアン演算には和、差、積があります。
https://ja.wikipedia.org/wiki/ブーリアン演算


DIM ST(2),C(3)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(ST,X,Y)
      IF DISTANCE(ST)<0 THEN
         CALL VEC3(C,1,0,0)
      ELSE
         CALL VEC3(C,0,0,0)
      END IF
      CALL SETCOLOR(C)
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  FUNCTION CIRCLE(P(),B(),R)
DIM PP(2)
MAT PP=P-B
LET CIRCLE=LENGTH(PP)-R
END FUNCTION

EXTERNAL  FUNCTION DISTANCE(P())
DIM P1(2),P2(2)
CALL VEC2(P1,-.25,0)
CALL VEC2(P2,.25,0)
LET A=CIRCLE(P,P1,.6)
LET B=CIRCLE(P,P2,.4)
!LET DISTANCE=A
!LET DISTANCE=B
LET DISTANCE=MIN(A,B)   ! ブーリアン演算 和
!LET DISTANCE=MAX(A,B)  ! ブーリアン演算 積
!LET DISTANCE=MAX(-A,B) ! ブーリアン演算 差
!LET DISTANCE=MAX(A,-B) ! ブーリアン演算 差
!LET DISTANCE=MIN(MAX(-A,B),MAX(A,-B)) ! ブーリアン演算
END FUNCTION

EXTERNAL  SUB SETCOLOR(COL())
SET COLOR COLORINDEX(CLAMP(COL(1),0,1),CLAMP(COL(2),0,1),CLAMP(COL(3),0,1))
END SUB

EXTERNAL  SUB VEC2(A(),X,Y)
LET A(1)=X
LET A(2)=Y
END SUB

EXTERNAL  SUB VEC3(A(),X,Y,Z)
LET A(1)=X
LET A(2)=Y
LET A(3)=Z
END SUB

EXTERNAL  FUNCTION CLAMP(X,A,B)
LET CLAMP=MIN(B,MAX(X,A))
END FUNCTION

EXTERNAL  FUNCTION LENGTH(A())
LET LENGTH=SQR(A(1)^2+A(2)^2)
END FUNCTION
 

フラグメントシェーダー 再帰呼び出し tree

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時28分42秒
返信・引用
  フラグメントシェーダー 再帰呼び出し tree

https://gam0022.net/blog/2017/03/02/raymarching-fold/
距離関数TREEはループとabs関数で再帰呼び出しの代わりをしています。

※シェーダー言語GLSLは再帰呼び出しができないのでループで代用します。


DIM ST(2),C(3)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(ST,X,Y)
      IF DISTANCE(ST)>0 THEN
         CALL VEC3(C,0,0,0)
      ELSE
         CALL VEC3(C,1,0,0)
      END IF
      CALL SETCOLOR(C(1),C(2),C(3))
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  FUNCTION BOX(P(),B()) !距離関数 BOX
DIM D(2),DD(2)
FOR I=1 TO 2
   LET D(I)=ABS(P(I))-B(I)
NEXT I
FOR I=1 TO 2
   LET DD(I)=MAX(D(I),0)
NEXT I
LET BOX=LENGTH(DD)+MIN(MAX(D(1),D(2)),0)
END FUNCTION

EXTERNAL  FUNCTION TREE(P()) !距離関数 TREE
DIM Q(2),SIZE(2)
LET SCALE=.8
CALL VEC2(SIZE,.1,.4)
LET D=BOX(P,SIZE)
MAT Q=P
FOR I=1 TO 7
   LET Q(1)=ABS(P(1))
   LET Q(2)=Q(2)-SIZE(2)
   CALL ROTATE(Q,-RAD(45))
   LET D=MIN(D,BOX(P,SIZE))
   MAT P=Q
   MAT SIZE=SCALE*SIZE
NEXT I
LET TREE=D
END FUNCTION

EXTERNAL  FUNCTION DISTANCE(P())
LET P(2)=P(2)+.5
LET DISTANCE=TREE(P)
END FUNCTION

EXTERNAL  SUB SETCOLOR(R,G,B)
SET COLOR COLORINDEX(CLAMP(R,0,1),CLAMP(G,0,1),CLAMP(B,0,1))
END SUB

EXTERNAL  SUB ROTATE(A(),ANGLE)
DIM M(2,2)
LET M(1,1)=COS(ANGLE)
LET M(1,2)=-SIN(ANGLE)
LET M(2,1)=SIN(ANGLE)
LET M(2,2)=COS(ANGLE)
MAT A=A*M
END SUB

EXTERNAL  FUNCTION STEP(A,X)
IF X<A THEN LET STEP=0 ELSE LET STEP=1
END FUNCTION

EXTERNAL  SUB VEC2(A(),X,Y)
LET A(1)=X
LET A(2)=Y
END SUB

EXTERNAL  SUB VEC3(A(),X,Y,Z)
LET A(1)=X
LET A(2)=Y
LET A(3)=Z
END SUB

EXTERNAL  FUNCTION CLAMP(X,A,B)
LET CLAMP=MIN(B,MAX(X,A))
END FUNCTION

EXTERNAL  FUNCTION LENGTH(A())
LET LENGTH=SQR(A(1)^2+A(2)^2)
END FUNCTION
 

フラグメントシェーダー 折り畳み abs関数

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時27分36秒
返信・引用
  フラグメントシェーダー 折り畳み abs関数

abs関数を使うと万華鏡のような効果をうみます。
(-1,-1)~(1,1)の範囲においてx=abs(x)とするとy軸対称になり
更にy=abs(y)とするとx軸対称になり万華鏡のようになります。

https://docs.google.com/presentation/d/12RrqyAkFanKmfL96ZHvhDCozE-_rKFPlU1YVwej4_bc/htmlpresent
※現在時刻(TIME)を使用しています。

DIM C(3),P(2),T(2)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
LET TI=INT(TIME)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(P,X,Y)
      FOR K=0 TO 4
         FOR I=1 TO 2
            LET P(I)=ABS(1.5*P(I))-1
         NEXT I
         CALL ROTATE(P,TI*K)
      NEXT K
      FOR I=1 TO 2
         LET T(I)=1-SMOOTHSTEP(.01,.02,ABS(P(I)))
         LET C(I)=MIX(P(I),1,T(1)+T(2))
      NEXT I
      CALL SETCOLOR(C(1),C(2),.5)
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  SUB ROTATE(A(),T)
DIM M(2,2)
LET S=SIN(T)
LET C=COS(T)
LET M(1,1)=C
LET M(1,2)=S
LET M(2,1)=-S
LET M(2,2)=C
MAT A=A*M
END SUB

EXTERNAL  SUB SETCOLOR(R,G,B)
SET COLOR COLORINDEX(CLAMP(R,0,1),CLAMP(G,0,1),CLAMP(B,0,1))
END SUB

EXTERNAL  SUB VEC2(A(),X,Y)
LET A(1)=X
LET A(2)=Y
END SUB

EXTERNAL  SUB VEC3(A(),X,Y,Z)
LET A(1)=X
LET A(2)=Y
LET A(3)=Z
END SUB

EXTERNAL  FUNCTION CLAMP(X,A,B)
LET CLAMP=MIN(B,MAX(X,A))
END FUNCTION

EXTERNAL  FUNCTION MIX(X,Y,A)
LET MIX=X*(1-A)+Y*A
END FUNCTION

EXTERNAL  FUNCTION SMOOTHSTEP(A,B,X)
LET T=CLAMP((X-A)/(B-A),0,1)
LET SMOOTHSTEP=T*T*(3-2*T)
END FUNCTION

EXTERNAL  FUNCTION SIGN(X)
IF X<0 THEN
   LET SIGN=-1
ELSEIF X>0 THEN
   LET SIGN=1
ELSE
   LET SIGN=0
END IF
END FUNCTION
 

フラグメントシェーダー 増殖 mod関数

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時26分30秒
返信・引用
  フラグメントシェーダー 増殖 mod関数

距離関数にmod関数を使うと増殖させることができます。

DIM ST(2),C(3)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(ST,X,Y)
      IF DISTANCE(ST)<0 THEN
         CALL VEC3(C,1,0,1)
      ELSE
         CALL VEC3(C,0,0,0)
      END IF
      CALL SETCOLOR(C(1),C(2),C(3))
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  FUNCTION CIRCLE(P(),R)
LET CIRCLE=LENGTH(P)-R
END FUNCTION

EXTERNAL  SUB REP(P(),C())
FOR I=1 TO 2
   LET P(I)=MOD(P(I)+.5*C(I),C(I))-.5*C(I)
NEXT I
END SUB

EXTERNAL  FUNCTION DISTANCE(P())
DIM T(2)
CALL VEC2(T,.25,.25)
CALL REP(P,T) !ここを注釈にすると1個だけになる。
LET DISTANCE=CIRCLE(P,.1)
END FUNCTION

EXTERNAL  SUB SETCOLOR(R,G,B)
SET COLOR COLORINDEX(CLAMP(R,0,1),CLAMP(G,0,1),CLAMP(B,0,1))
END SUB

EXTERNAL  SUB VEC2(A(),X,Y)
LET A(1)=X
LET A(2)=Y
END SUB

EXTERNAL  SUB VEC3(A(),X,Y,Z)
LET A(1)=X
LET A(2)=Y
LET A(3)=Z
END SUB

EXTERNAL  FUNCTION CLAMP(X,A,B)
LET CLAMP=MIN(B,MAX(X,A))
END FUNCTION

EXTERNAL  FUNCTION LENGTH(A())
LET LENGTH=SQR(A(1)^2+A(2)^2)
END FUNCTION
 

フラグメントシェーダー ハート

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時25分29秒
返信・引用
  フラグメントシェーダー ハート

ハート型を定義しています。

DIM P(2),COL(3)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(P,X,Y)
      !! MAT P=2*P
      IF DISTANCE(P)<0 THEN
         CALL VEC3(COL,1,0,0)
      ELSE
         CALL VEC3(COL,0,0,0)
      END IF
      CALL SETCOLOR(COL(1),COL(2),COL(3))
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  FUNCTION HEART(P()) !距離関数 ハート
DIM P1(2),P2(2),P3(2)
CALL VEC2(P1,.25,.75)
CALL VEC2(P2,0,1)
LET P(1)=ABS(P(1))
IF P(2)+P(1)>1 THEN
   MAT P1=P-P1
   LET HEART=SQR(DOT(P1,P1))-SQR(2)/4
ELSE
   MAT P2=P-P2
   FOR I=1 TO 2
      LET P3(I)=P(I)-.5*MAX(P(1)+P(2),0)
   NEXT I
   LET HEART=SQR(MIN(DOT(P2,P2),DOT(P3,P3)))*SIGN(P(1)-P(2))
END IF
END FUNCTION

EXTERNAL  FUNCTION DISTANCE(P())
DIM P1(2)
CALL VEC2(P1,0,-.5)
MAT P=P-P1
LET DISTANCE=HEART(P)
!!LET DISTANCE=(P(1)^2+P(2)^2-1)^3-P(1)^2*P(2)^3
END FUNCTION

EXTERNAL  SUB SETCOLOR(R,G,B)
SET COLOR COLORINDEX(CLAMP(R,0,1),CLAMP(G,0,1),CLAMP(B,0,1))
END SUB

EXTERNAL  SUB VEC2(A(),X,Y)
LET A(1)=X
LET A(2)=Y
END SUB

EXTERNAL  FUNCTION SIGN(X)
IF X<0 THEN
   LET SIGN=-1
ELSEIF X>0 THEN
   LET SIGN=1
ELSE
   LET SIGN=0
END IF
END FUNCTION

EXTERNAL  SUB VEC3(A(),X,Y,Z)
LET A(1)=X
LET A(2)=Y
LET A(3)=Z
END SUB

EXTERNAL  FUNCTION CLAMP(X,A,B)
LET CLAMP=MIN(B,MAX(X,A))
END FUNCTION
 

フラグメントシェーダー 星

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時24分54秒
返信・引用
  フラグメントシェーダー 星

星型を定義しています。


DIM ST(2),C(3)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(ST,X,Y)
      IF DISTANCE(ST)>0 THEN
         CALL VEC3(C,0,0,0)
      ELSE
         CALL VEC3(C,1,0,0)
      END IF
      CALL SETCOLOR(C(1),C(2),C(3))
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  FUNCTION STAR5(P(),R,RF) !距離関数 星型
DIM K1(2),K2(2),BA(2),P1(2),PP(2)
CALL VEC2(K1,.809016994375,-.587785252292)
CALL VEC2(K2,-K1(1),K1(2))
LET P(1)=ABS(P(1))
MAT PP=P
FOR I=1 TO 2
   LET P(I)=P(I)-2*MAX(DOT(K1,PP),0)*K1(I)
NEXT I
MAT PP=P
FOR I=1 TO 2
   LET P(I)=P(I)-2*MAX(DOT(K2,PP),0)*K2(I)
NEXT I
LET P(1)=ABS(P(1))
LET P(2)=P(2)-R
CALL VEC2(K2,-K1(2),K1(1))
MAT BA=RF*K2
CALL VEC2(K2,0,1)
MAT BA=BA-K2
LET H=CLAMP(DOT(P,BA)/DOT(BA,BA),0,R)
MAT P1=H*BA
MAT P1=P-P1
LET STAR5=LENGTH(P1)*SIGN(P(2)*BA(1)-P(1)*BA(2))
END FUNCTION

EXTERNAL  FUNCTION DISTANCE(P())
LET DISTANCE=STAR5(P,.5,.381966011250105)
END FUNCTION

EXTERNAL  SUB SETCOLOR(R,G,B)
SET COLOR COLORINDEX(CLAMP(R,0,1),CLAMP(G,0,1),CLAMP(B,0,1))
END SUB

EXTERNAL  SUB VEC2(A(),X,Y)
LET A(1)=X
LET A(2)=Y
END SUB

EXTERNAL  SUB VEC3(A(),X,Y,Z)
LET A(1)=X
LET A(2)=Y
LET A(3)=Z
END SUB

EXTERNAL  FUNCTION CLAMP(X,A,B)
LET CLAMP=MIN(B,MAX(X,A))
END FUNCTION

EXTERNAL  FUNCTION LENGTH(A())
LET LENGTH=SQR(A(1)^2+A(2)^2)
END FUNCTION

EXTERNAL  FUNCTION SIGN(X)
IF X<0 THEN
   LET SIGN=-1
ELSEIF X>0 THEN
   LET SIGN=1
ELSE
   LET SIGN=0
END IF
END FUNCTION
 

フラグメントシェーダー 正多角形

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時24分7秒
返信・引用
  フラグメントシェーダー 正多角形

正多角形を定義しています。
http://math.sakura.ne.jp/?action=common_download_main&upload_id=2052


DIM P(2),COL(3)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(P,X,Y)
      IF DISTANCE(P)<0 THEN
         CALL VEC3(COL,1,0,0)
      ELSE
         CALL VEC3(COL,0,0,0)
      END IF
      CALL SETCOLOR(COL(1),COL(2),COL(3))
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  FUNCTION POLYGON(P(),N)
FOR K=1 TO INT(N/2)
   LET S=S+ABS(P(1)*SIN(2*K*PI/N)-ABS(P(2))*COS(2*K*PI/N))
   LET E=E+SIN(2*K*PI/N)
NEXT K
LET S=S+.5*ABS(P(2)*SIN(N*PI/2))-.5*P(1)*COT(INT(N/2)*PI/N)
LET E=E-.5*COT(INT(N/2)*PI/N)
LET POLYGON=S-E
END FUNCTION

EXTERNAL  FUNCTION DISTANCE(P())
CALL ROTATE(P,PI/6)
LET DISTANCE=POLYGON(P,3)

!CALL ROTATE(P,PI/4)
!LET DISTANCE=POLYGON(P,4)

!CALL ROTATE(P,PI/5*1.5)
!LET DISTANCE=POLYGON(P,5)

!LET DISTANCE=POLYGON(P,6)

!CALL ROTATE(P,PI/14)
!LET DISTANCE=POLYGON(P,7)

!CALL ROTATE(P,PI/8)
!LET DISTANCE=POLYGON(P,8)
END FUNCTION

EXTERNAL  SUB ROTATE(A(),T)
DIM M(2,2)
LET S=SIN(T)
LET C=COS(T)
LET M(1,1)=C
LET M(1,2)=S
LET M(2,1)=-S
LET M(2,2)=C
MAT A=A*M
END SUB

EXTERNAL  SUB SETCOLOR(R,G,B)
SET COLOR COLORINDEX(CLAMP(R,0,1),CLAMP(G,0,1),CLAMP(B,0,1))
END SUB

EXTERNAL  SUB VEC2(A(),X,Y)
LET A(1)=X
LET A(2)=Y
END SUB

EXTERNAL  SUB VEC3(A(),X,Y,Z)
LET A(1)=X
LET A(2)=Y
LET A(3)=Z
END SUB

EXTERNAL  FUNCTION CLAMP(X,A,B)
LET CLAMP=MIN(B,MAX(X,A))
END FUNCTION
 

フラグメントシェーダー 五角形

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時23分15秒
返信・引用
  フラグメントシェーダー 五角形

五角形を定義しています。


DIM ST(2),C(3)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(ST,X,Y)
      IF DISTANCE(ST)>0 THEN
         CALL VEC3(C,0,0,0)
      ELSE
         CALL VEC3(C,1,0,0)
      END IF
      CALL SETCOLOR(C(1),C(2),C(3))
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  FUNCTION PENTAGON(P(),R) !距離関数 五角形
DIM K(3),P1(2),P2(2),PP(2)
CALL VEC3(K,.809016994,.587785252,.726542528)
CALL VEC2(P1,-K(1),K(2))
CALL VEC2(P2,K(1),K(2))
LET P(1)=ABS(P(1))
MAT PP=P
FOR I=1 TO 2
   LET P(I)=P(I)-2*MIN(DOT(P1,PP),0)*P1(I)
NEXT I
MAT PP=P
FOR I=1 TO 2
   LET P(I)=P(I)-2*MIN(DOT(P2,PP),0)*P2(I)
NEXT I
CALL VEC2(P1,CLAMP(P(1),-R*K(3),R*K(3)),R)
MAT P=P-P1
LET PENTAGON=LENGTH(P)*SIGN(P(2))
END FUNCTION

EXTERNAL  FUNCTION DISTANCE(P())
LET DISTANCE=PENTAGON(P,.5)
!!LET DISTANCE=ABS(PENTAGON(P,.5))-.1 !くり抜き
!!LET DISTANCE=2*ABS(P(2))+ABS(SQR(3)*P(1)-P(2))+ABS(SQR(3)*P(1)+P(2))-2*SQR(3) !6角形
END FUNCTION

EXTERNAL  SUB SETCOLOR(R,G,B)
SET COLOR COLORINDEX(CLAMP(R,0,1),CLAMP(G,0,1),CLAMP(B,0,1))
END SUB

EXTERNAL  SUB VEC2(A(),X,Y)
LET A(1)=X
LET A(2)=Y
END SUB

EXTERNAL  SUB VEC3(A(),X,Y,Z)
LET A(1)=X
LET A(2)=Y
LET A(3)=Z
END SUB

EXTERNAL  FUNCTION CLAMP(X,A,B)
LET CLAMP=MIN(B,MAX(X,A))
END FUNCTION

EXTERNAL  FUNCTION LENGTH(A())
LET LENGTH=SQR(A(1)^2+A(2)^2)
END FUNCTION

EXTERNAL  FUNCTION SIGN(X)
IF X<0 THEN
   LET SIGN=-1
ELSEIF X>0 THEN
   LET SIGN=1
ELSE
   LET SIGN=0
END IF
END FUNCTION
 

フラグメントシェーダー 三角形

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時22分28秒
返信・引用
  フラグメントシェーダー 三角形

三角形を定義しています。

DIM P(2),COL(3),PP(2)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(P,X,Y)
      IF DISTANCE(P)<0 THEN
         CALL VEC3(COL,1,0,0)
      ELSE
         CALL VEC3(COL,0,0,0)
      END IF
      CALL SETCOLOR(COL(1),COL(2),COL(3))
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  FUNCTION TRIANGLE(P(),P0(),P1(),P2()) !距離関数 三角形
DIM E0(2),E1(2),E2(2),V0(2),V1(2),V2(2),PQ0(2),PQ1(2),PQ2(2),D(2)
MAT E0=P1-P0
MAT E1=P2-P1
MAT E2=P0-P2
MAT V0=P-P0
MAT V1=P-P1
MAT V2=P-P2
FOR I=1 TO 2
   LET PQ0(I)=V0(I)-E0(I)*CLAMP(DOT(V0,E0)/DOT(E0,E0),0,1)
   LET PQ1(I)=V1(I)-E1(I)*CLAMP(DOT(V1,E1)/DOT(E1,E1),0,1)
   LET PQ2(I)=V2(I)-E2(I)*CLAMP(DOT(V2,E2)/DOT(E2,E2),0,1)
NEXT I
LET S=SIGN(E0(1)*E2(2)-E0(2)*E2(1))
CALL VEC2(V0,DOT(PQ0,PQ0),S*(V0(1)*E0(2)-V0(2)*E0(1)))
CALL VEC2(V1,DOT(PQ1,PQ1),S*(V1(1)*E1(2)-V1(2)*E1(1)))
CALL VEC2(V2,DOT(PQ2,PQ2),S*(V2(1)*E2(2)-V2(2)*E2(1)))
FOR I=1 TO 2
   LET D(I)=MIN(MIN(V0(I),V1(I)),V2(I))
NEXT I
LET TRIANGLE=-SQR(D(1))*SIGN(D(2))
END FUNCTION

EXTERNAL  FUNCTION DISTANCE(P())
DIM P0(2),P1(2),P2(2)
CALL VEC2(P0,-.9,.1)
CALL VEC2(P1,.3,.5)
CALL VEC2(P2,-.5,-.8)
LET DISTANCE=TRIANGLE(P,P0,P1,P2)
!!LET DISTANCE=TRIANGLE(P,P0,P1,P2)-.1 ! 角を丸くする
!!LET DISTANCE=ABS(4*P(1)+(2*SQR(3))*ABS(P(2))-1)+2*SQR(3)*ABS(P(2))-3
END FUNCTION

EXTERNAL  SUB SETCOLOR(R,G,B)
SET COLOR COLORINDEX(CLAMP(R,0,1),CLAMP(G,0,1),CLAMP(B,0,1))
END SUB

EXTERNAL  SUB VEC2(A(),X,Y)
LET A(1)=X
LET A(2)=Y
END SUB

EXTERNAL  SUB VEC3(A(),X,Y,Z)
LET A(1)=X
LET A(2)=Y
LET A(3)=Z
END SUB

EXTERNAL  FUNCTION CLAMP(X,A,B)
LET CLAMP=MIN(B,MAX(X,A))
END FUNCTION

EXTERNAL  FUNCTION SIGN(X)
IF X<0 THEN
   LET SIGN=-1
ELSEIF X>0 THEN
   LET SIGN=1
ELSE
   LET SIGN=0
END IF
END FUNCTION
 

フラグメントシェーダー 箱

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時21分48秒
返信・引用
  フラグメントシェーダー 箱

距離関数BOXに四角形(箱)を定義しています。
どうしてこれで四角形(箱)が描けるのかを考えてみるのもいいかもしれません。

また距離関数ROUNDBOXは角に丸みがある四角形(箱)です。
距離関数BOXとの違いはRを引いているかだけです。

DIM ST(2),C(3)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(ST,X,Y)
      IF DISTANCE(ST)>0 THEN
         CALL VEC3(C,0,0,0)
      ELSE
         CALL VEC3(C,1,0,0)
      END IF
      CALL SETCOLOR(C(1),C(2),C(3))
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  FUNCTION BOX(P(),B()) !距離関数 箱
DIM D(2),DD(2)
FOR I=1 TO 2
   LET D(I)=ABS(P(I))-B(I)
NEXT I
FOR I=1 TO 2
   LET DD(I)=MAX(D(I),0)
NEXT I
LET BOX=LENGTH(DD)+MIN(MAX(D(1),D(2)),0)
END FUNCTION

EXTERNAL  FUNCTION ROUNDBOX(P(),B(),R) !距離関数
LET ROUNDBOX=BOX(P,B)-R
END FUNCTION

EXTERNAL  FUNCTION DISTANCE(P())
DIM B(2)
CALL VEC2(B,.3,.3)
LET DISTANCE=BOX(P,B)
!!LET DISTANCE=ROUNDBOX(P,B,.1)
!!LET DISTANCE=ABS((P(1)-B(1))+(P(2)-B(2)))+ABS((P(1)-B(1))-(P(2)-B(2)))-1
END FUNCTION

EXTERNAL  SUB SETCOLOR(R,G,B)
SET COLOR COLORINDEX(CLAMP(R,0,1),CLAMP(G,0,1),CLAMP(B,0,1))
END SUB

EXTERNAL  SUB VEC2(A(),X,Y)
LET A(1)=X
LET A(2)=Y
END SUB

EXTERNAL  SUB VEC3(A(),X,Y,Z)
LET A(1)=X
LET A(2)=Y
LET A(3)=Z
END SUB

EXTERNAL  FUNCTION CLAMP(X,A,B)
LET CLAMP=MIN(B,MAX(X,A))
END FUNCTION

EXTERNAL  FUNCTION LENGTH(A())
LET LENGTH=SQR(A(1)^2+A(2)^2)
END FUNCTION
 

フラグメントシェーダー 円

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時20分54秒
返信・引用
  フラグメントシェーダー 円

ここでは距離関数に円を定義しています。
これは円の方程式 x^2+y^2=r^2 から SQR(x^2+y^2)-r=0となります。
距離関数が正か負で塗分けています。

では距離関数をどう定義したらいいのかについては
具体的な距離関数がネット上に公開されています。

https://iquilezles.org/www/articles/distfunctions2d/distfunctions2d.htm
https://www.shadertoy.com/playlist/MXdSRf

3次元の距離関数も公開されています。
https://www.iquilezles.org/www/articles/distfunctions/distfunctions.htm


DIM ST(2),C(3)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE)
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE)
      CALL VEC2(ST,X,Y)
      IF DISTANCE(ST)>0 THEN
         CALL VEC3(C,0,0,0)
      ELSE
         CALL VEC3(C,1,0,0)
      END IF
      CALL SETCOLOR(C(1),C(2),C(3))
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  FUNCTION CIRCLE(P(),R) !距離関数 円
LET CIRCLE=LENGTH(P)-R
!!LET CIRCLE=NLENGTH(P,4)-R
END FUNCTION

EXTERNAL  FUNCTION ELLIPSE(P(),A,B,N) !距離関数 楕円
LET ELLIPSE=ABS(P(1)/A)^N+ABS(P(2)/B)^N-1
END FUNCTION

EXTERNAL  FUNCTION DISTANCE(P())
LET DISTANCE=CIRCLE(P,.5)
!LET DISTANCE=ELLIPSE(P,.7,.4,2)
END FUNCTION

EXTERNAL  SUB SETCOLOR(R,G,B)
SET COLOR COLORINDEX(CLAMP(R,0,1),CLAMP(G,0,1),CLAMP(B,0,1))
END SUB

EXTERNAL  SUB VEC2(A(),X,Y)
LET A(1)=X
LET A(2)=Y
END SUB

EXTERNAL  SUB VEC3(A(),X,Y,Z)
LET A(1)=X
LET A(2)=Y
LET A(3)=Z
END SUB

EXTERNAL  FUNCTION CLAMP(X,A,B)
LET CLAMP=MIN(B,MAX(X,A))
END FUNCTION

EXTERNAL  FUNCTION LENGTH(A())
LET LENGTH=SQR(A(1)^2+A(2)^2)
END FUNCTION

EXTERNAL  FUNCTION NLENGTH(P(),N)
LET NLENGTH=(ABS(P(1))^N+ABS(P(2))^N)^(1/N)
END FUNCTION
 

フラグメントシェーダー グラフ

 投稿者:しばっち  投稿日:2022年 2月11日(金)14時20分11秒
返信・引用
  フラグメントシェーダー グラフ

https://6317.teacup.com/basic/bbs/4548
https://6317.teacup.com/basic/bbs/4549
https://6317.teacup.com/basic/bbs/4550
https://6317.teacup.com/basic/bbs/4551
https://6317.teacup.com/basic/bbs/4552



OpenGLのシェーダー言語GLSLをもとにしています。
フラグメントシェーダーでは各ドットに対して色を決めて描画していくもので
使用する描画コマンドはPLOT POINTSのみです。PLOT LINESやDRAW文は使いません。

各R,G,Bを0~1までの数値で表したフルカラーを使います。
WINDOWは(-1,-1)-(1,1)として原点(0,0)を中心とした正方形領域を想定しています。


フラグメントシェーダーにはMOD関数による繰り返しやABS関数による折り畳み
ブーリアン演算などがあります。


フラグメントシェーダーでは基本的に距離関数というものを定義して描画していきます。
その距離関数の正体は陰関数(f(x,y)=0)です。y=f(x)もy-f(x)=0とすれば陰関数になります。
陰関数なのでf(x,y)<0ならその内側になり、f(x,y)>0なら外側になります。

ここではy=3*sin(2*x)^3+sin(x)^2-2*sin(3*x)のグラフを描画しています。


なお、シェーダー言語GLSL(WebGL)ではGPUをごりごりに使用したアニメーション表示です。(GPUで暖がとれる!?)
見ているだけでも楽しいのでぜひアクセスしてみてください。

https://glslsandbox.com/
https://glslsandbox.com/e#78473.0
https://glslsandbox.com/e#78672.0
https://glslsandbox.com/e#78754.0
https://glslsandbox.com/e#78747.0
https://glslsandbox.com/e#78666.0
https://glslsandbox.com/e#78263.0
https://glslsandbox.com/e#78701.0
https://glslsandbox.com/e#78390.0


DIM P(2),COL(3),CC(3),B(3)
ASK BITMAP SIZE XSIZE,YSIZE
SET COLOR MODE "NATIVE"
SET POINT STYLE 1
SET WINDOW -XSIZE/MIN(XSIZE,YSIZE),XSIZE/MIN(XSIZE,YSIZE),-YSIZE/MIN(XSIZE,YSIZE),YSIZE/MIN(XSIZE,YSIZE) !! SET WINDOW -1,1,-1,1
CALL VEC3(CC,1,.4,.1)
FOR YY=0 TO YSIZE-1
   FOR XX=0 TO XSIZE-1
      LET X=(XX*2-XSIZE)/MIN(XSIZE,YSIZE) ! -1~1
      LET Y=(YY*2-YSIZE)/MIN(XSIZE,YSIZE) ! -1~1
      CALL VEC2(P,X,Y)
      MAT P=5*P
      IF ABS(MOD(P(1),1))<.02 OR ABS(MOD(P(2),1))<.02 THEN CALL VEC3(B,.2,.2,.2) ELSE CALL VEC3(B,0,0,0)
      IF ABS(P(1))<.03 OR ABS(P(2))<.03 THEN CALL VEC3(B,.5,.5,.5)
      LET L=ABS(DISTANCE(P)-.6)
      IF L<>0 THEN
         LET F=.1/L
      ELSE
         LET F=1
      END IF
      MAT COL=F*CC
      MAT COL=COL+B
      CALL SETCOLOR(COL(1),COL(2),COL(3))
      PLOT POINTS:X,Y
   NEXT XX
NEXT YY
END

EXTERNAL  FUNCTION DISTANCE(P()) !距離関数 y=3*sin(2*x)^3+sin(x)^2-2*sin(3*x)
LET DISTANCE=3*SIN(2*P(1))^3+SIN(P(1))^2-2*SIN(3*P(1))-P(2)
END FUNCTION

EXTERNAL  SUB SETCOLOR(R,G,B)
SET COLOR COLORINDEX(CLAMP(R,0,1),CLAMP(G,0,1),CLAMP(B,0,1))
END SUB

EXTERNAL  SUB VEC2(A(),X,Y)
LET A(1)=X
LET A(2)=Y
END SUB

EXTERNAL  SUB VEC3(A(),X,Y,Z)
LET A(1)=X
LET A(2)=Y
LET A(3)=Z
END SUB

EXTERNAL  FUNCTION CLAMP(X,A,B)
LET CLAMP=MIN(B,MAX(X,A))
END FUNCTION
 

オイラーのゼータ関数

 投稿者:たろさ  投稿日:2022年 1月16日(日)21時36分47秒
返信・引用
  ! 6n+1,6n-1 篩 Ver.2
!5億までの素数を配列DIMに登録する。
!オイラーのゼータ関数
OPTION ARITHMETIC NATIVE       !2進モード
LET t0=TIME
LET Pk=26355867
LET k=5E8    !1E7  9999991 素数 (664579th)
LET k1=IP(k/6)
LET k2=IP(SQR(k1))

DIM A5(k1),A7(k1)

MAT A5=ZER
MAT A7=ZER

LET P1=5
LET C1=1

DO
   FOR n=1 TO k1
      LET P6=P1*n+C1
      IF P6 > k1 THEN EXIT FOR
      LET A5(P6)=1
   NEXT n
   FOR n=1 TO k1
      LET P6=P1*n-C1
      IF P6 > k1 THEN EXIT FOR
      LET A7(P6)=1
   NEXT n

   LET P1=P1+6
   LET C1=C1+1
   IF P1 >k1 THEN EXIT DO
LOOP


LET P1=7
LET C1=1
DO
   FOR n=1 TO k1
      LET P6=P1*n+C1
      IF P6 > k1 THEN EXIT FOR
      LET A7(P6)=1
   NEXT n

   FOR n=1 TO k1
      LET P6=P1*n-C1
      IF P6 > k1 THEN EXIT FOR
      LET A5(P6)=1
   NEXT n

   LET P1=P1+6
   LET C1=C1+1
   IF P1 >k1 THEN EXIT DO
LOOP
DIM B(Pk)       !1E8 99999989 素数 (5761455th)
LET B(1)=2      !2E8 199999991 素数 (11078937th)
LET B(2)=3
LET C=3
FOR n=1 TO k1
   LET P5=6*n-1
   LET P7=6*n+1
   IF A5(n)=0 THEN
      LET B(C)=P5
      LET C=C+1
   END IF
   IF A7(n)=0 THEN
      LET B(C)=P7
      LET C=C+1
   END IF
NEXT n

LET Q=1
FOR n=1 TO Pk
   LET P=B(n)^2
   LET Q=Q*(P/(P-1))
NEXT n

PRINT Q
PRINT PI^2/6

LET TM=TIME-t0
PRINT TM;"秒"

END

動作確認 Microsoft surface pro6  メモリーDDR3 8GB   Win11 Pro 64-bit
Paract BASIC Ver. 2.1.3.1 (2022.01.09)
64ビット版Lazarus3.2.2

計算結果
1.644934066247
1.644934066848
9.474000000002 秒

5億まで では?
15億までは確認しました。1.644934066

とある八雲の科学解説 『素数から生まれる円周率』
https://www.youtube.com/watch?v=aq_-cIv92t8
 

オイラーの列の和を求める方法

 投稿者:たろさ  投稿日:2022年 1月16日(日)20時51分21秒
返信・引用 編集済
  解析教程〈上〉 単行本 ? 1997/10/1
E. ハイラー  (著), G. ワナー  (著)

1.5.3 サイン関数とオイラー積 P74~75
(5.23)   1/pi^2 +1/4pi^2 +1/9pi^2 +1/16pi^2 +1/25pi^2+ … =1/6


FOR n=1 TO 1E9
   LET Z=Z+1/(PI^2*n^2)
NEXT n
PRINT Z
PRINT 1/6
END


計算結果
.1666666656523
.1666666666667


動作確認 Microsoft surface pro6  メモリーDDR3 8GB   Win11 Pro 64-bit
Paract BASIC Ver. 2.1.3.1 (2022.01.09)
64ビット版Lazarus3.2.2


e^(pi*i)=-1  ;https://keisan.casio.jp/calculator

PRINT EXP(PI*COMPLEX(0,1))
計算結果
(-1  1.224646799147E-16)
(-1  1.22460635382238E-16)


https://keisan.casio.jp/exec/system/1161228679
zeta(2)=pi^2/6
zeta(100)=pi^100*189196075638244250590454866138987443745405683066133872266771392408622790830394495422/9815205420757514710108178059369553458327392260750404049930407987933582359080767225644716670683512153512547802166033089160919189453125
https://6317.teacup.com/basic/bbs/4115
 

Re: 白石先生へ

 投稿者:RCカー  投稿日:2022年 1月11日(火)08時51分37秒
返信・引用
  SHIRAISHI Kazuoさんへのお返事です。
ありがとうございます。
> RCカーさんへのお返事です。
>
> 十進BASICの行列計算(MAT文)は厳密解を求めるアルゴリズムになっています。
> 十進BASICのMAT文はサイズが大きい行列の計算には不向きなので,行列の特徴に応じたアルゴリズムに組み直した方がベターです。大きい行列の計算には反復計算で収束させるようなアルゴリズムがあるので,そのような例が見出される本を選べばいいと思います(お勧めは思いつきません)。
>
 

Re: 白石先生へ

 投稿者:SHIRAISHI Kazuo  投稿日:2022年 1月11日(火)08時10分0秒
返信・引用
  > No.4999[元記事へ]

RCカーさんへのお返事です。

十進BASICの行列計算(MAT文)は厳密解を求めるアルゴリズムになっています。
十進BASICのMAT文はサイズが大きい行列の計算には不向きなので,行列の特徴に応じたアルゴリズムに組み直した方がベターです。大きい行列の計算には反復計算で収束させるようなアルゴリズムがあるので,そのような例が見出される本を選べばいいと思います(お勧めは思いつきません)。

http://hp.vector.co.jp/authors/VA008683/

 

白石先生へ

 投稿者:RCカー  投稿日:2022年 1月10日(月)15時29分30秒
返信・引用
  十進ベーシックでの、行列計算について、参考になるものはありませんか。
最小二乗法で、観測方程式と、正規方程式と、連立方程式の解き方、χ二乗検定について、
教えてほしいです。
観測値の数が、多すぎると、余剰観測といって、解が一つに定まりません。
少ないと、解けません。余剰観測値がゼロになるとそれはそれで、困るようです。
西先生の本と原田先生の本です。
最小二乗法は、測量以外にも応用が利くようですとのこと。
 

Re: 十進BASICのバグだとは思いますが...

 投稿者:SHIRAISHI Kazuo  投稿日:2022年 1月 2日(日)08時13分51秒
返信・引用 編集済
  10進モードの十進演算ルーチンは指数部をかなり大きく取っています。
べき乗の計算精度がほぼ確実と考えられる範囲をもとにMAXNUMとEPS(0)を定めています。
Full BASICの十進演算は,数値式と数値変数は別物です。
MAXNUMとEPS(0)は,数値変数の精度からきまる定数です。
十進BASICだと,数値式の途中でMAXNUMを超えるような計算も可能です。

10 LET x=MAXNUM
20 PRINT (MAXNUM*1000)/10000
30 END
数値変数にEPS(0)より小さい正の数を代入すると,変数の値は0になります。
十進BASICでは,数値式の中でそれより小さい数も扱えます。
10 LET x=EPS(0)
20 PRINT (x/100)*1000
30 END

10進1000桁モードのEPS(0)は1E-1017です。
10 OPTION ARITHMETIC DECIMAL_HIGH
20 PRINT USING "#.##^^^^^^":EPS(0)
30 END

ただし,MAXNUMより大きい数やEPS(0)より小さい数を数値変数に保持できるのは,十進BASICのバグかも知れません。








> WINDOWS版十進BASIC ver. 7.8.5.6で下記のプログラムを
> 2進モード、複素数モードで実行するとI=171でエラーとなり実行が止まります。
> 10進モードだとI=70で止まります。
> 1000桁モードではI=453で止まります。
>
> LET A=1
> FOR I=2 TO 500
>    LET A=A*I
>    PRINT I;A
> NEXT I
> END
>
> 次に下記のようにプログラムを変更して2進モード、複素数モードで実行すると
> I=170以降は同じ値が表示されます。
>
> ですが10進モードだとそのまま計算し続けます。
> しかもMAXNUM(1E99)をはるかに超えています。
>
> WHEN EXCEPTION IN文がエラーによる停止を抑止したため、MAXNUMを超えて計算しています。
> また、MAXNUM値を超えた数値が正しく判定できています。
>
> これは明らかに十進BASICのバグだとは思いますが、2進モード(複素数モード)での限界(1.79769313486232E308)を超えていますので
> 指数表記ではありますが、正しく計算できているようなので超巨大数の計算に利用できそうです。(巨大数の計算には有理数モードがありますが...)
>
> 但し、MAXNUM値を大幅に超えていますので何らかの不具合を起こす可能性があるかもしれません。
> (実際にどこまで計算できるのかは試しておりません。2000の階乗までは問題ない!?)
>
> 2進モード(複素数モード)ではIEEE754による数値表現で限界値があるのはわかりますが
> 10進モードでの数値表現(どんなものかよくわかりませんが)では、まだまだ計算可能ということだと思います。
>
>
> 1000桁モードではI=453でMAXNUMを超えています。
> I=819で小数点以下2000桁を超えています。
> I=820から指数表示に変わりますがそのまま計算し続けてます。
>
> LET A=1
> FOR I=2 TO 2000
>
>    WHEN EXCEPTION IN
>       LET A=A*I
>    USE
>    END WHEN
>
>    PRINT I;
>    IF MAXNUM<A THEN PRINT "TRUE"; ELSE PRINT "FALSE";
>    PRINT A
> NEXT I
> END
>
>        2進モード(複素数モード)での実行(一部のみ)
>
>  165 FALSE 5.42391066613159E295
>  166 FALSE 9.00369170577843E297
>  167 FALSE 1.503616514865E300
>  168 FALSE 2.5260757449732E302
>  169 FALSE 4.2690680090047E304
>  170 FALSE 7.25741561530799E306
>  171 FALSE 7.25741561530799E306
>  172 FALSE 7.25741561530799E306
>  173 FALSE 7.25741561530799E306
>
>        10進モードでの実行(一部のみ)
>
>  65 FALSE 8.24765059208253E90
>  66 FALSE 5.44344939077447E92
>  67 FALSE 3.64711109181889E94
>  68 FALSE 2.48003554243685E96
>  69 FALSE 1.71122452428143E98
>  70 TRUE 1.197857166997E100
>  71 TRUE 8.50478588567871E101
>  72 TRUE 6.12344583768867E103
>  73 TRUE 4.47011546151273E105
>  74 TRUE 3.30788544151942E107
>  75 TRUE 2.48091408113956E109
>
>  1997 TRUE 4.151569143494E5725
>  1998 TRUE 8.29483514870102E5728
>  1999 TRUE 1.65813754622533E5732
>  2000 TRUE 3.31627509245067E5735
>
>
> ちなみに下記のようにした場合は
> 2進モード、複素数モードで実行するとI=178以降で0になります。
>
> 10進モードではI=70以降で0になります。
> (I=70以降も0にはならずにずっと指数表示が続くかと思ったのですが...)
>
> 1000桁モードではI=456で小数点以下2000桁を超えていますが
> I=457で突然0になります。これは何らかの処理により0が代入されたものと思いますが
> この挙動はいささか不自然にも思えます。
>
> LET A=1
> FOR I=2 TO 500
>    LET A=A/I
>    PRINT I;A
> NEXT I
> END
>
 

十進BASICのバグだとは思いますが...

 投稿者:しばっち  投稿日:2022年 1月 1日(土)19時05分44秒
返信・引用
  WINDOWS版十進BASIC ver. 7.8.5.6で下記のプログラムを
2進モード、複素数モードで実行するとI=171でエラーとなり実行が止まります。
10進モードだとI=70で止まります。
1000桁モードではI=453で止まります。

LET A=1
FOR I=2 TO 500
   LET A=A*I
   PRINT I;A
NEXT I
END

次に下記のようにプログラムを変更して2進モード、複素数モードで実行すると
I=170以降は同じ値が表示されます。

ですが10進モードだとそのまま計算し続けます。
しかもMAXNUM(1E99)をはるかに超えています。

WHEN EXCEPTION IN文がエラーによる停止を抑止したため、MAXNUMを超えて計算しています。
また、MAXNUM値を超えた数値が正しく判定できています。

これは明らかに十進BASICのバグだとは思いますが、2進モード(複素数モード)での限界(1.79769313486232E308)を超えていますので
指数表記ではありますが、正しく計算できているようなので超巨大数の計算に利用できそうです。(巨大数の計算には有理数モードがありますが...)

但し、MAXNUM値を大幅に超えていますので何らかの不具合を起こす可能性があるかもしれません。
(実際にどこまで計算できるのかは試しておりません。2000の階乗までは問題ない!?)

2進モード(複素数モード)ではIEEE754による数値表現で限界値があるのはわかりますが
10進モードでの数値表現(どんなものかよくわかりませんが)では、まだまだ計算可能ということだと思います。


1000桁モードではI=453でMAXNUMを超えています。
I=819で小数点以下2000桁を超えています。
I=820から指数表示に変わりますがそのまま計算し続けてます。

LET A=1
FOR I=2 TO 2000

   WHEN EXCEPTION IN
      LET A=A*I
   USE
   END WHEN

   PRINT I;
   IF MAXNUM<A THEN PRINT "TRUE"; ELSE PRINT "FALSE";
   PRINT A
NEXT I
END

       2進モード(複素数モード)での実行(一部のみ)

165 FALSE 5.42391066613159E295
166 FALSE 9.00369170577843E297
167 FALSE 1.503616514865E300
168 FALSE 2.5260757449732E302
169 FALSE 4.2690680090047E304
170 FALSE 7.25741561530799E306
171 FALSE 7.25741561530799E306
172 FALSE 7.25741561530799E306
173 FALSE 7.25741561530799E306

       10進モードでの実行(一部のみ)

65 FALSE 8.24765059208253E90
66 FALSE 5.44344939077447E92
67 FALSE 3.64711109181889E94
68 FALSE 2.48003554243685E96
69 FALSE 1.71122452428143E98
70 TRUE 1.197857166997E100
71 TRUE 8.50478588567871E101
72 TRUE 6.12344583768867E103
73 TRUE 4.47011546151273E105
74 TRUE 3.30788544151942E107
75 TRUE 2.48091408113956E109

1997 TRUE 4.151569143494E5725
1998 TRUE 8.29483514870102E5728
1999 TRUE 1.65813754622533E5732
2000 TRUE 3.31627509245067E5735


ちなみに下記のようにした場合は
2進モード、複素数モードで実行するとI=178以降で0になります。

10進モードではI=70以降で0になります。
(I=70以降も0にはならずにずっと指数表示が続くかと思ったのですが...)

1000桁モードではI=456で小数点以下2000桁を超えていますが
I=457で突然0になります。これは何らかの処理により0が代入されたものと思いますが
この挙動はいささか不自然にも思えます。

LET A=1
FOR I=2 TO 500
   LET A=A/I
   PRINT I;A
NEXT I
END
 

日付時刻の国際標準表記

 投稿者:nagram  投稿日:2021年12月31日(金)17時23分41秒
返信・引用
  今年もあとわずかなので、日付時刻に関するプログラムです。
日付や時刻の表記法は様々な様式がありますが、ISO 8601 (日本では JIS X 0301) で規格化されています。
規格に沿った日付時刻を表示する関数を作りました。
皆さま、良いお年をお迎えください。

DECLARE EXTERNAL FUNCTION datetime1$, datetimesec1$, datetime2$, datetimesec2$
!
PRINT "基本形式"
PRINT datetime1$      ! "YYYYMMDDThhmmss+0900"
PRINT datetimesec1$   ! "YYYYMMDDThhmmss.ss+0900"
PRINT "拡張形式"
PRINT datetime2$      ! "YYYY-MM-DDThh:mm:ss+09:00"
PRINT datetimesec2$   ! "YYYY-MM-DDThh:mm:ss.ss+09:00"
END
!
EXTERNAL FUNCTION datetime1$  ! 日付時刻国際標準表記(基本形式)
LET d0$=DATE$
LET t$=TIME$
LET d$=DATE$
IF d$<>d0$ THEN LET t$="00:00:00"  ! 日付更新考慮
LET datetime1$=d$&"T"&t$(1:2)&t$(4:5)&t$(7:8)&"+0900"
END FUNCTION
!
EXTERNAL FUNCTION datetimesec1$  ! 日付時刻国際標準表記(基本形式,小数秒)
LET t0=TIME
LET d$=DATE$
LET t=TIME
IF t<t0 THEN LET d$=DATE$  ! 日付更新考慮
LET h=INT(t/3600)
LET hm$=RIGHT$("0"&STR$(h),2)&RIGHT$("0"&STR$(INT((t-3600*h)/60)),2)
LET t$=USING$("%%.##",MOD(t,60))
! LET t1=MOD(t,60)  ! 秒の小数部末尾の"0"を記述しない
! IF t1>=10 THEN LET t$=STR$(t1) ELSE LET t$="0"&STR$(t1)
! IF t1<1 THEN LET t$="0"&t$
LET datetimesec1$=d$&"T"&hm$&t$&"+0900"
END FUNCTION
!
EXTERNAL FUNCTION datetime2$  ! 日付時刻国際標準表記(拡張形式)
LET d0$=DATE$
LET t$=TIME$
LET d$=DATE$
IF d$<>d0$ THEN LET t$="00:00:00"  ! 日付更新考慮
LET datetime2$=d$(1:4)&"-"&d$(5:6)&"-"&d$(7:8)&"T"&t$&"+09:00"
END FUNCTION
!
EXTERNAL FUNCTION datetimesec2$  ! 日付時刻国際標準表記(拡張形式,小数秒)
LET t0=TIME
LET d$=DATE$
LET t=TIME
IF t<t0 THEN LET d$=DATE$  ! 日付更新考慮
LET h=INT(t/3600)
LET hm$=RIGHT$("0"&STR$(h),2)&":"&RIGHT$("0"&STR$(INT((t-3600*h)/60)),2)
LET t$=USING$("%%.##",MOD(t,60))
! LET t1=MOD(t,60)  ! 秒の小数部末尾の"0"を記述しない
! IF t1>=10 THEN LET t$=STR$(t1) ELSE LET t$="0"&STR$(t1)
! IF t1<1 THEN LET t$="0"&t$
LET datetimesec2$=d$(1:4)&"-"&d$(5:6)&"-"&d$(7:8)&"T"&hm$&":"&t$&"+09:00"
END FUNCTION
 

白石先生へ

 投稿者:RCカー  投稿日:2021年12月19日(日)22時01分36秒
返信・引用
  有難う御座います  

Re: 御助言お願いします の2

 投稿者:SHIRAISHI Kazuo  投稿日:2021年12月19日(日)20時22分33秒
返信・引用
  > No.4992[元記事へ]

RCカーさんへのお返事です。

> 十進ベーシックはウィンドウズ11対応予定はありますか?
>
おそらく動作するものと思います。
不具合があればお知らせください。
 

Re: 御助言お願いします

 投稿者:SHIRAISHI Kazuo  投稿日:2021年12月19日(日)20時21分39秒
返信・引用
  > No.4991[元記事へ]

RCカーさんへのお返事です。

> 構造化ベーシックについて今現在購入できる。
> お勧め図書の御推薦をお待ちしています。
>
中古限定となりますが,入手可能です。

数学とコンピュータ1



数学とコンピュータ2
 森北出版 新数学入門シリーズ
 

御助言お願いします の2

 投稿者:RCカー  投稿日:2021年12月19日(日)05時21分53秒
返信・引用
  十進ベーシックはウィンドウズ11対応予定はありますか?
 

御助言お願いします

 投稿者:RCカー  投稿日:2021年12月19日(日)05時20分15秒
返信・引用
  構造化ベーシックについて今現在購入できる。
お勧め図書の御推薦をお待ちしています。
 

DOT言語変換

 投稿者:しばっち  投稿日:2021年12月12日(日)10時01分48秒
返信・引用
  隣接行列からDOT言語に変換します。
隣接行列からグラフを作成するのはなかなか大変な作業ですが、このDOT言語からグラフ画像を作成してくれます。

DOT言語自体は無向グラフ、有向グラフに対応していますが、下記プログラムは有向グラフ版です。("digraph"を"graph"に、"->"を"--"に変更するだけですが...)
実行するとdotファイル(DOT言語)を書き出します。

このdotファイルをメモ帳で開いてもたいした意味はありませんので
下記URLよりダウンロードしてください。私はwindows版のzip版(win32)をダウンロードしました。

http://www.graphviz.org/


dot.exeがあるフォルダ(bin)でコマンドプロンプトから

dot sample.dot -Tpng -o sample.png

と打ち込むとグラフ画像(pngファイル)が出来上がります。

https://ja.wikipedia.org/wiki/DOT言語
http://cbrc3.cbrc.jp/~tominaga/translations/graphviz/dotguide.pdf
https://qiita.com/rubytomato@github/items/51779135bc4b77c8c20d


DIM NODE$(26)
MAT READ NODE$
DATA A,B,C,D,E,F,G,H,I,J,K,L,M,N,O,P,Q,R,S,T,U,V,W,X,Y,Z
LET F$="D:\tool\Graphviz\bin\sample" !ドライブ、パスを指定してください。(ファイル名はsample.dotになります)
LET FORMAT$="png"  ! jpg gif svg ps pdf...
OPEN #1:NAME F$&".dot"
ERASE #1
PRINT #1:"digraph sample {"
!!PRINT #1:"  graph[rankdir=LR];" ! TB BT LR RL
!!PRINT #1:"  node [shape=box];"  ! box polygon ellipse circle point egg triangle diamond trapezium parallellogram house pentagon hexagon doublecircle...
RESTORE 10
READ SIZE
FOR J=1 TO SIZE
   FOR I=1 TO SIZE
      READ S$
      IF S$<>"INF" THEN
         PRINT #1:"    ";NODE$(J);" -> ";NODE$(I);" [fontcolor=red label=";CHR$(34);S$;CHR$(34);"];"
      END IF
   NEXT I
NEXT J
PRINT #1:"}"
CLOSE #1
!!EXECUTE "D:\tool\Graphviz\bin\dot.exe" WITH(F$&".dot","-T"&FORMAT$,"-o",F$&"."&FORMAT$)
!!GLOAD F$&"."&FORMAT$

10 DATA 10 ! 有向グラフ
   DATA INF,  1,  4,INF,  5,INF,  2,INF,  2,INF
   DATA INF,INF,  3,INF,INF,  4,INF,INF,INF,INF
   DATA   3,INF,INF,  2,INF,INF,  1,INF,INF,  1
   DATA INF,INF,  2,INF,INF,INF,INF,  3,INF,INF
   DATA INF,  2,INF,INF,INF,  2,INF,INF,INF,  1
   DATA   1,INF,  4,INF,  1,INF,  5,INF,INF,INF
   DATA INF,  3,INF,INF,INF,  6,INF,  5,INF,INF
   DATA   1,INF,INF,  3,INF,INF,INF,INF,  5,INF
   DATA INF,  2,INF,  5,INF,  2,INF,  3,INF,  2
   DATA   1,INF,  4,  3,INF,INF,  5,INF,INF,INF

20 DATA 8 ! 無向グラフ(対角成分に対して対称。A(I,J)とA(J,I)が同値)
   DATA INF,  1,  1,INF,INF,INF,  1,  1
   DATA   1,INF,  1,INF,  1,INF,INF,INF
   DATA   1,  1,INF,  1,INF,  1,  1,INF
   DATA INF,INF,  1,INF,  1,INF,  1,  1
   DATA INF,  1,INF,  1,INF,  1,INF,INF
   DATA INF,INF,  1,INF,  1,INF,  1,INF
   DATA   1,INF,  1,  1,INF,  1,INF,  1
   DATA   1,INF,INF,  1,INF,INF,  1,INF

30 DATA 10 ! 有向グラフ
   DATA INF,  2,INF,  1,INF,  4,  6,INF,  7,INF
   DATA INF,INF,INF, -2,  1,INF,INF,  2,INF,  1
   DATA  -3,INF,INF,INF,  1,INF, -1,INF,  2,INF
   DATA INF,  2,INF,INF,INF,INF,  1,INF,INF,  3
   DATA INF,INF, -1,INF,INF,  2,INF,  3,  4,INF
   DATA   4,INF,INF,INF,  2,INF,INF,  1,INF,  5
   DATA INF,  1,INF,  2,INF,  1,INF,  5,INF,INF
   DATA   2,INF,INF,  4,INF,  3,INF,INF,  2,  4
   DATA INF,INF,  4,INF,  4,INF, -2,INF,INF,INF
   DATA   1,INF,  3,INF,  6,INF,  1,INF,INF,INF
END

 

ダイクストラ法

 投稿者:しばっち  投稿日:2021年12月12日(日)09時58分39秒
返信・引用
  https://ja.wikipedia.org/wiki/ダイクストラ法

ダイクストラ法によるグラフ探索です。
このプログラムは移植版です。


DIM NODE$(26)
MAT READ NODE$
DATA A,B,C,D,E,F,G,H,I,J,K,L,M,N,O,P,Q,R,S,T,U,V,W,X,Y,Z
RESTORE 10
READ SIZE
DIM VISITED(SIZE),DISTANCE(SIZE),INDEX(SIZE),A(SIZE,SIZE)
LET TRUE=1
LET FALSE=0
LET INF=100000000
MAT DISTANCE=(INF)*CON
FOR I=1 TO SIZE
   FOR J=1 TO SIZE
      READ S$
      IF S$<>"INF" THEN LET A(I,J)=VAL(S$) ELSE LET A(I,J)=INF
   NEXT J
NEXT I
LET START=3
LET DISTANCE(START)=0
LET NEXTINDEX=START
DO
   LET I=NEXTINDEX
   LET VISITED(I)=TRUE
   LET LMIN=INF
   FOR J=1 TO SIZE
      IF VISITED(J)=FALSE THEN
         IF A(I,J)<>INF AND DISTANCE(J)>DISTANCE(I)+A(I,J) THEN
            LET DISTANCE(J)=DISTANCE(I)+A(I,J)
            LET INDEX(J)=I
         END IF
         IF DISTANCE(J)<LMIN THEN
            LET LMIN=DISTANCE(J)
            LET NEXTINDEX=J
         END IF
      END IF
   NEXT J
LOOP WHILE LMIN<INF
FOR I=1 TO SIZE
   IF I<>START AND VISITED(I)=TRUE THEN
      LET K=I
      PRINT NODE$(K);
      DO
         LET K=INDEX(K)
         PRINT " ← ";NODE$(K);
      LOOP UNTIL K=START
      PRINT "  DISTANCE";DISTANCE(I)
   END IF
NEXT I

10 DATA 10 ! 有向グラフ
   DATA INF,  1,  4,INF,  5,INF,  2,INF,  2,INF
   DATA INF,INF,  3,INF,INF,  4,INF,INF,INF,INF
   DATA   3,INF,INF,  2,INF,INF,  1,INF,INF,  1
   DATA INF,INF,  2,INF,INF,INF,INF,  3,INF,INF
   DATA INF,  2,INF,INF,INF,  2,INF,INF,INF,  1
   DATA   1,INF,  4,INF,  1,INF,  5,INF,INF,INF
   DATA INF,  3,INF,INF,INF,  6,INF,  5,INF,INF
   DATA   1,INF,INF,  3,INF,INF,INF,INF,  5,INF
   DATA INF,  2,INF,  5,INF,  2,INF,  3,INF,  2
   DATA   1,INF,  4,  3,INF,INF,  5,INF,INF,INF

20 DATA 8 ! 無向グラフ(対角成分に対して対称。A(I,J)とA(J,I)が同値)
   DATA INF,  1,  1,INF,INF,INF,  1,  1
   DATA   1,INF,  1,INF,  1,INF,INF,INF
   DATA   1,  1,INF,  1,INF,  1,  1,INF
   DATA INF,INF,  1,INF,  1,INF,  1,  1
   DATA INF,  1,INF,  1,INF,  1,INF,INF
   DATA INF,INF,  1,INF,  1,INF,  1,INF
   DATA   1,INF,  1,  1,INF,  1,INF,  1
   DATA   1,INF,INF,  1,INF,INF,  1,INF

30 DATA 10 ! 有向グラフ
   DATA INF,  2,INF,  1,INF,  4,  6,INF,  7,INF
   DATA INF,INF,INF, -2,  1,INF,INF,  2,INF,  1
   DATA  -3,INF,INF,INF,  1,INF, -1,INF,  2,INF
   DATA INF,  2,INF,INF,INF,INF,  1,INF,INF,  3
   DATA INF,INF, -1,INF,INF,  2,INF,  3,  4,INF
   DATA   4,INF,INF,INF,  2,INF,INF,  1,INF,  5
   DATA INF,  1,INF,  2,INF,  1,INF,  5,INF,INF
   DATA   2,INF,INF,  4,INF,  3,INF,INF,  2,  4
   DATA INF,INF,  4,INF,  4,INF, -2,INF,INF,INF
   DATA   1,INF,  3,INF,  6,INF,  1,INF,INF,INF

40 DATA 7
   DATA INF,  1,INF,  3,INF,  5,INF
   DATA   2,INF,  3,INF,INF,INF,  2
   DATA INF,  2,INF,  2,  4,INF,INF
   DATA INF,INF,  4,INF,INF,INF,  5
   DATA   3,INF,  3,  3,INF,  2,INF
   DATA INF,  1,INF,  4,  2,INF,INF
   DATA   5,INF,  2,INF,  3,  6,INF
END
 

深さ優先グラフ探索

 投稿者:しばっち  投稿日:2021年12月12日(日)09時57分19秒
返信・引用
  深さ優先の再帰呼び出しによるグラフ探索です。


PUBLIC NUMERIC A(20,20),DISTANCE(20),VISITED(20),INDEX(100),SIZE,INF,TRUE,FALSE,START
PUBLIC NUMERIC MININDEX(100),MINDISTANCE,MAXINDEX(100),MAXDISTANCE
PUBLIC STRING NODE$(26)
MAT READ NODE$
DATA A,B,C,D,E,F,G,H,I,J,K,L,M,N,O,P,Q,R,S,T,U,V,W,X,Y,Z
LET INF=100000000
LET MINDISTANCE=INF
LET MAXDISTANCE=-INF
LET TRUE=1
LET FALSE=0
RESTORE 10
READ SIZE
FOR I=1 TO SIZE ! 隣接行列読み込み
   FOR J=1 TO SIZE
      READ S$
      IF S$="INF" THEN
         LET A(I,J)=INF
      ELSE
         LET A(I,J)=VAL(S$)
      END IF
   NEXT J
NEXT I
LET START=8
LET GOAL=2
IF START<1 OR GOAL<1 OR START>SIZE OR GOAL>SIZE OR START=GOAL THEN
   PRINT "ERROR !!"
   STOP
END IF
CALL VISIT(START,GOAL)
PRINT REPEAT$("-",60)
LET K=START
DO
   PRINT NODE$(K);" → ";
   LET K=MININDEX(K)
LOOP UNTIL K=GOAL
PRINT NODE$(GOAL);"  MIN DISTANCE";MINDISTANCE
!LET K=START
!DO
!   PRINT NODE$(K);" → ";
!   LET K=MAXINDEX(K)
!LOOP UNTIL K=GOAL
!PRINT NODE$(GOAL);"  MAX DISTANCE";MAXDISTANCE

10 DATA 10 ! 有向グラフ
   DATA INF,  1,  4,INF,  5,INF,  2,INF,  2,INF
   DATA INF,INF,  3,INF,INF,  4,INF,INF,INF,INF
   DATA   3,INF,INF,  2,INF,INF,  1,INF,INF,  1
   DATA INF,INF,  2,INF,INF,INF,INF,  3,INF,INF
   DATA INF,  2,INF,INF,INF,  2,INF,INF,INF,  1
   DATA   1,INF,  4,INF,  1,INF,  5,INF,INF,INF
   DATA INF,  3,INF,INF,INF,  6,INF,  5,INF,INF
   DATA   1,INF,INF,  3,INF,INF,INF,INF,  5,INF
   DATA INF,  2,INF,  5,INF,  2,INF,  3,INF,  2
   DATA   1,INF,  4,  3,INF,INF,  5,INF,INF,INF

20 DATA 8 ! 無向グラフ(対角成分に対して対称。A(I,J)とA(J,I)が同値)
   DATA INF,  1,  1,INF,INF,INF,  1,  1
   DATA   1,INF,  1,INF,  1,INF,INF,INF
   DATA   1,  1,INF,  1,INF,  1,  1,INF
   DATA INF,INF,  1,INF,  1,INF,  1,  1
   DATA INF,  1,INF,  1,INF,  1,INF,INF
   DATA INF,INF,  1,INF,  1,INF,  1,INF
   DATA   1,INF,  1,  1,INF,  1,INF,  1
   DATA   1,INF,INF,  1,INF,INF,  1,INF

30 DATA 10 ! 有向グラフ
   DATA INF,  2,INF,  1,INF,  4,  6,INF,  7,INF
   DATA INF,INF,INF, -2,  1,INF,INF,  2,INF,  1
   DATA  -3,INF,INF,INF,  1,INF, -1,INF,  2,INF
   DATA INF,  2,INF,INF,INF,INF,  1,INF,INF,  3
   DATA INF,INF, -1,INF,INF,  2,INF,  3,  4,INF
   DATA   4,INF,INF,INF,  2,INF,INF,  1,INF,  5
   DATA INF,  1,INF,  2,INF,  1,INF,  5,INF,INF
   DATA   2,INF,INF,  4,INF,  3,INF,INF,  2,  4
   DATA INF,INF,  4,INF,  4,INF, -2,INF,INF,INF
   DATA   1,INF,  3,INF,  6,INF,  1,INF,INF,INF
END

EXTERNAL  SUB VISIT(I,GOAL)
   IF I=GOAL THEN
      LET K=START
      DO
         PRINT NODE$(K);" → ";
         LET K=INDEX(K)
      LOOP UNTIL K=GOAL
      PRINT NODE$(GOAL);"  DISTANCE";DISTANCE(GOAL)
      IF MINDISTANCE>DISTANCE(GOAL) THEN
         MAT MININDEX=INDEX
         LET MINDISTANCE=DISTANCE(GOAL)
      END IF
      !      IF MAXDISTANCE<DISTANCE(GOAL) THEN
      !         MAT MAXINDEX=INDEX
      !         LET MAXDISTANCE=DISTANCE(GOAL)
      !      END IF
   ELSE
   !IF MINDISTANCE>DISTANCE(I) THEN
      FOR J=1 TO SIZE
         IF I<>J AND A(I,J)<>INF AND VISITED(J)=FALSE THEN
            LET VISITED(I)=TRUE
            LET DISTANCE(J)=DISTANCE(I)+A(I,J)
            LET INDEX(I)=J
            CALL VISIT(J,GOAL)
            LET VISITED(I)=FALSE
            LET INDEX(I)=0
            LET DISTANCE(J)=0
         END IF
      NEXT J
      !END IF
   END IF
END SUB
 

ベルマンフォード法

 投稿者:永野護  投稿日:2021年11月30日(火)12時24分19秒
返信・引用
  nagram様よりご提示いただきましたプログラムをヒントにdo  while   loopを加えてみました。
隣接行列の値も変えています。これでちゃんとできているみたいですが。

DIM  dist(7), node$(7),path$(7)
MAT READ node$
DATA  A,B,C,D,E,F,G
MAT path$=node$(1)&NUL$

DIM  c(7,7)
MAT c=1000*CON
LET inf=1000
LET dist(1)=0
LET dist(2)=inf
LET dist(3)=inf
LET dist(4)=inf
LET dist(5)=inf
LET dist(6)=inf
LET dist(7)=inf
LET c(1,1)=0
LET c(2,2)=0
LET c(3,3)=0
LET c(4,4)=0
LET c(5,5)=0
LET c(6,6)=0
LET c(7,7)=0
LET c(1,2)=30
LET c(2,1)=30
LET c(1,3)=20
LET c(3,1)=20
LET c(1,7)=1
LET c(7,1)=1
LET c(2,3)=2
LET c(3,2)=2
LET c(2,4)=6
LET c(4,2)=6
LET c(2,5)=1
LET c(5,2)=1
LET c(2,7)=6
LET c(7,2)=6
LET c(3,4)=2
LET c(4,3)=2
LET c(3,5)=3
LET c(5,3)=3
LET c(3,6)=1
LET c(6,3)=1
LET c(4,5)=20
LET c(5,4)=20
LET c(4,6)=30
LET c(6,4)=30
LET c(4,7)=100
LET c(7,4)=100
LET c(5,7)=2
LET c(7,5)=2
LET c(5,6)=5
LET c(6,5)=5
LET c(6,7)=10
LET c(7,6)=10

REM-----------------------------------------------------------------
LET w=0
DO  WHILE  w<=7
   FOR  t=1  TO 7
      FOR i=1  TO 7

         IF  dist(i)>dist(t)+c(t,i)  THEN
            LET dist(i)=dist(t)+c(t,i)
            LET path$(i)=path$(t)&node$(i)


         END IF

      NEXT  i

   NEXT  t
   LET w=w+1
LOOP
REM-----------------------------------
PRINT  "-------------------------------------------------------"

FOR   x=1  TO 7
   PRINT  "最短距離="; dist(x); ", 経路=";path$(x)
NEXT x

END





 

ベルマンフォード法

 投稿者:永野護  投稿日:2021年11月29日(月)12時31分35秒
返信・引用
  nagram様、お忙しいところ貴重なプログラムを作っていただきましたことに心より感謝
します。ご多忙中にもかかわらず丁寧なご指導をいただき、ありがとうございました。
貴重な体験を糧とし、日々精進してまいりたいと存じます。
末筆ではございますが、nagram様の一層のご活躍を心よりお祈り申し上げます。
敬具
 

Re: ベルマンフォード法

 投稿者:nagram  投稿日:2021年11月29日(月)11時09分35秒
返信・引用
  永野護さんへのお返事です。

> ベルマンフォード法は無向グラフにも適用できるのでしょうか。先ごろ作っていただきました
> プログラムで対応可能でしょうか。度々質問してすみませんがご迷惑でなかったらご教授
> お願いしたいのですが。何卒よろしくお願いします。本当に身勝手なお願いで申し訳ないです。
>

私もグラフ理論に詳しいわけではないですが、ベルマンフォード法は有向グラフにしか使えないようですね。
無向グラフの解法も調べましたが、いい方法は見つけられませんでした。

永野護さんが提示したプログラムも有向グラフにしか適用できず、「辺の方向は頂点の順序により決定する」という条件が付きます。
c(p,q) を p<q のときのみ調査しています。(A→B,D→G は調査するが、E→C はスルー)
提示されたグラフのデータがその条件を満たしていたので解けましたが、c(5,2)=6 といったデータがあると上手く求めることができません。

ネットでベルマンフォード法のサンプルプログラムを見ましたが、[FOR i=1 TO 7] の部分を foreach文を利用しているものが多いですね。
foreach文は十進BASICにはない構文なので代替する方法を考えましたが、私の実力不足でうまくできませんでした。
下記のプログラムは、永野護さんが提示したプログラムを探査する向きを変えて2回実行しています。当然、計算量も2倍になります。
(データの一部を変えて、頂点の順序とは逆方向の辺を作っています。)
この方法が、(非負の重みの)すべての有向グラフに対して有効かは分かりません。


各頂点をA~Gとし、有向グラフの重みを次のようにし、Aを始点として各頂点への最短経路を求める。

A→B=20, A→C=50
B→C=150, B→D=6
C→E=4
D→C=3, D→E=8, D→G=20
E→F=40, E→G=30
F→G=80

LET n=7
DIM  dist(n), node$(n),path$(n)
MAT READ node$
DATA  A,B,C,D,E,F,G
MAT path$=node$(1)&NUL$

DIM  c(n,n)
MAT c=1000*CON

MAT dist=1000*CON
LET dist(1)=0

LET c(1,2)=20
LET c(1,3)=50
LET c(2,3)=150  ! c(2,3)=1
LET c(2,4)=6
LET c(4,3)=3  ! c(3,4)=100
LET c(3,5)=4  ! c(3,5)=50
LET c(4,5)=8
LET c(4,7)=20
LET c(5,7)=30
LET c(5,6)=40
LET c(6,7)=80
REM-----------------------------------------------------------------
FOR  t=1  TO n-1
   FOR i=1  TO  n
      IF     dist(i)>dist(t)+c(t,i)  THEN
         LET dist(i)=dist(t)+c(t,i)
         LET path$(i)=path$(t)&node$(i)
      END IF
   NEXT  i
NEXT  t
FOR  t=1  TO n-1
   FOR i=n  TO  1 STEP -1
      IF     dist(i)>dist(t)+c(t,i)  THEN
         LET dist(i)=dist(t)+c(t,i)
         LET path$(i)=path$(t)&node$(i)
      END IF
   NEXT  i
NEXT  t

REM   ------------------------------

FOR   x=1  TO  n
   PRINT  "最短距離="; dist(x); ", 経路=";path$(x)
NEXT x

END
 

ベルマンフォード法

 投稿者:永野護  投稿日:2021年11月27日(土)18時31分44秒
返信・引用
  ベルマンフォード法は無向グラフにも適用できるのでしょうか。先ごろ作っていただきました
プログラムで対応可能でしょうか。度々質問してすみませんがご迷惑でなかったらご教授
お願いしたいのですが。何卒よろしくお願いします。本当に身勝手なお願いで申し訳ないです。
 

(無題)

 投稿者:永野護  投稿日:2021年11月25日(木)11時09分21秒
返信・引用
  nagramさん、丁寧な解説ありがとうございました。
助かりました。
 

Re: ベルマンフォード法

 投稿者:nagram  投稿日:2021年11月25日(木)10時18分33秒
返信・引用
  > No.4981[元記事へ]

永野護さんへのお返事です。

各頂点をA~Gとし、有向グラフの重みを次のようにし、Aを始点として各頂点への最短経路を求める問題ですね。

A→B=20, A→C=50
B→C=1, B→D=6
C→D=100, C→E=50
D→E=8, D→G=20
E→F=40, E→G=30
F→G=80


DIM  dist(7), node$(7),path$(7)
MAT READ node$
DATA  A,B,C,D,E,F,G
MAT path$=node$(1)&NUL$

DIM  c(7,7)
MAT c=1000*CON

LET dist(1)=0
LET dist(2)=1000
LET dist(3)=1000
LET dist(4)=1000
LET dist(5)=1000
LET dist(6)=1000
LET dist(7)=1000
LET c(1,1)=0
LET c(2,2)=0
LET c(3,3)=0
LET c(4,4)=0
LET c(5,5)=0
LET c(6,6)=0
LET c(7,7)=0
LET c(1,2)=20
LET c(2,1)=20
LET c(1,3)=50
LET c(3,1)=50
LET c(2,3)=1
LET c(3,2)=1
LET c(2,4)=6
LET c(4,2)=6
LET c(3,4)=100
LET c(4,3)=100
LET c(3,5)=50
LET c(5,3)=50
LET c(4,5)=8
LET c(5,4)=8
LET c(4,7)=20
LET c(7,4)=20
LET c(5,7)=30
LET c(7,5)=30
LET c(5,6)=40
LET c(6,5)=40
LET c(6,7)=80
LET c(7,6)=80
REM-----------------------------------------------------------------
FOR  t=1  TO 6
   FOR i=1  TO  7
      IF     dist(i)>dist(t)+c(t,i)  THEN
         LET dist(i)=dist(t)+c(t,i)
         LET path$(i)=path$(t)&node$(MAX(t,i))
      END IF
   NEXT  i
NEXT  t
REM   ------------------------------

FOR   x=1  TO  7
   PRINT  "最短距離="; dist(x); ", 経路=";path$(x)
NEXT x

END
 

ベルマンフォード法

 投稿者:永野護  投稿日:2021年11月24日(水)10時40分51秒
返信・引用
  REM   ベルマンフォード法で最短経路を求めています。最短距離は何とか出てるみたいですが、途中経路の
REM   示し方がわかりません。ご多忙の折誠に恐縮ですが、よろしければ経路の示し方をご教授願えないでしょうか。
REM   コードは下記の通りです。7つの頂点を含むグラフです。

DIM  dist(7)

DIM  c(7,7)
FOR  y=1  TO  7
   FOR  z=1  TO  7
      LET c(y,z)=1000
   NEXT Z
NEXT Y

LET dist(1)=0
LET dist(2)=1000
LET dist(3)=1000
LET dist(4)=1000
LET dist(5)=1000
LET dist(6)=1000
LET dist(7)=1000
LET c(1,1)=0
LET c(2,2)=0
LET c(3,3)=0
LET c(4,4)=0
LET c(5,5)=0
LET c(6,6)=0
LET c(7,7)=0
LET c(1,2)=20
LET c(2,1)=20
LET c(1,3)=50
LET c(3,1)=50
LET c(2,3)=1
LET c(3,2)=1
LET c(2,4)=6
LET c(4,2)=6
LET c(3,4)=100
LET c(4,3)=100
LET c(3,5)=50
LET c(5,3)=50
LET c(4,5)=8
LET c(5,4)=8
LET c(4,7)=20
LET c(7,4)=20
LET c(5,7)=30
LET c(7,5)=30
LET c(5,6)=40
LET c(6,5)=40
LET c(6,7)=80
LET c(7,6)=80
REM-----------------------------------------------------------------
FOR  t=1  TO 6
   FOR i=1  TO  7
    IF     dist(i)>dist(t)+c(t,i)  THEN  LET dist(i)=dist(t)+c(t,i)
  NEXT  i
NEXT  t

REM   ------------------------------


FOR   x=1  TO  7
   PRINT  "最短距離="; dist(x);
NEXT x

END










 

1億~1兆までの素数の個数を数える計算時間 14分21.37秒

 投稿者:たろさ  投稿日:2021年11月 8日(月)10時54分52秒
返信・引用 編集済
  !30n+k篩  8スレッド     2021/07/10
!
!https://6317.teacup.com/basic/bbs/4360
!
!Paract BASIC 30n+k篩 Ver.12  500兆  5/9  (1E8) step
!
DECLARE STRUCTURE struct1: 1 OF NUMERIC
DECLARE STRUCTURE struct2: 1 OF NUMERIC(1536277)
DECLARE STRUCTURE struct4: 3 OF NUMERIC
DECLARE SHARED sha OF struct2
DECLARE MESSAGE mes2 OF struct4
DECLARE MESSAGE mes3 OF struct4
DECLARE MESSAGE mes4 OF struct4
DECLARE MESSAGE mes5 OF struct4
DECLARE MESSAGE mes6 OF struct4
DECLARE MESSAGE mes7 OF struct4
DECLARE MESSAGE mes8 OF struct4
DECLARE MESSAGE met2 OF struct1
DECLARE MESSAGE met3 OF struct1
DECLARE MESSAGE met4 OF struct1
DECLARE MESSAGE met5 OF struct1
DECLARE MESSAGE met6 OF struct1
DECLARE MESSAGE met7 OF struct1
DECLARE MESSAGE met8 OF struct1
PARACT PART1
OPTION ARITHMETIC NATIVE
LET t1=TIME
LET k=24494963
LET k3=1536277
DECLARE EXTERNAL SUB prime
CALL prime(k)
WAIT EVENT Ok5
LET S=1E8  !110E11  !pi(1E12),37607912018
LET E=1E10 !111E11  !pi(1E11),4118054813    (1E10)455052511
LET ST=1E8

START Part2
START PART3
START PART4
START PART5
START PART6
START PART7
START PART8
SEND TO mes2 FROM S,E,ST
SEND TO mes3 FROM S,E,ST
SEND TO mes4 FROM S,E,ST
SEND TO mes5 FROM S,E,ST
SEND TO mes6 FROM S,E,ST
SEND TO mes7 FROM S,E,ST
SEND TO mes8 FROM S,E,ST
LET TOTAL=5761455
DECLARE EXTERNAL FUNCTION cprime
!DECLARE NUMERIC X,Y,Z

FOR I=S TO E-ST STEP ST
   LET t0=TIME
   LET L=cprime(I,I+ST/8)
   RECEIVE FROM met2 TO X
   RECEIVE FROM met3 TO Y
   RECEIVE FROM met4 TO Z
   RECEIVE FROM met5 TO X1
   RECEIVE FROM met6 TO Y1
   RECEIVE FROM met7 TO Z1
   RECEIVE FROM met8 TO X2
   LET L=L+X+Y+Z+X1+Y1+Z1+X2
   LET TOTAL=TOTAL+L
   !IF MOD(I+ST,1E9)=0 THEN PRINT (I+ST)/1E9;TOTAL
   !PRINT TOTAL
   PRINT (I+ST)/1E8;TOTAL;L;
   LET TM=TIME-t0
   PRINT USING"###.###":TM;
   PRINT "秒"
NEXT I
LET TM=TIME-t1
PRINT USING"#####.##":TM;
PRINT "秒"
END PARACT

PARACT PART2
OPTION ARITHMETIC NATIVE
RECEIVE FROM mes2 TO S,E,ST
DECLARE EXTERNAL FUNCTION cprime
DECLARE NUMERIC L
FOR I=S TO E-ST STEP ST
   LET L=cprime(I+ST/8,I+ST/4)
   SEND TO met2 FROM L
NEXT I
END PARACT

PARACT PART3
OPTION ARITHMETIC NATIVE
RECEIVE FROM mes3 TO S,E,ST
DECLARE EXTERNAL FUNCTION cprime
DECLARE NUMERIC L
FOR I=S TO E-ST STEP ST
   LET L=cprime(I+ST/4,I+3*ST/8)
   SEND TO met3 FROM L
NEXT I
END PARACT

PARACT PART4
OPTION ARITHMETIC NATIVE
RECEIVE FROM mes4 TO S,E,ST
DECLARE EXTERNAL FUNCTION cprime
DECLARE NUMERIC L
FOR I=S TO E-ST STEP ST
   LET L=cprime(I+3*ST/8,I+ST/2)
   SEND TO met4 FROM L
NEXT I
END PARACT

PARACT PART5
OPTION ARITHMETIC NATIVE
RECEIVE FROM mes5 TO S,E,ST
DECLARE EXTERNAL FUNCTION cprime
DECLARE NUMERIC L
FOR I=S TO E-ST STEP ST
   LET L=cprime(I+ST/2,I+5*ST/8)
   SEND TO met5 FROM L
NEXT I
END PARACT

PARACT PART6
OPTION ARITHMETIC NATIVE
RECEIVE FROM mes6 TO S,E,ST
DECLARE EXTERNAL FUNCTION cprime
DECLARE NUMERIC L
FOR I=S TO E-ST STEP ST
   LET L=cprime(I+5*ST/8,I+3*ST/4)
   SEND TO met6 FROM L
NEXT I
END PARACT

PARACT PART7
OPTION ARITHMETIC NATIVE
RECEIVE FROM mes7 TO S,E,ST
DECLARE EXTERNAL FUNCTION cprime
DECLARE NUMERIC L
FOR I=S TO E-ST STEP ST
   LET L=cprime(I+3*ST/4,I+7*ST/8)
   SEND TO met7 FROM L
NEXT I
END PARACT

PARACT PART8
OPTION ARITHMETIC NATIVE
RECEIVE FROM mes8 TO S,E,ST
DECLARE EXTERNAL FUNCTION cprime
DECLARE NUMERIC L
FOR I=S TO E-ST STEP ST
   LET L=cprime(I+7*ST/8,I+ST)
   SEND TO met8 FROM L
NEXT I
END PARACT

EXTERNAL FUNCTION cprime(k4,k6)
OPTION ARITHMETIC NATIVE
DECLARE NUMERIC G(1536277) !素数
GET FROM sha TO G

DIM B(8)
DATA 1,7,11,13,17,19,23,29
MAT READ B
LET Q=30
LET U=IP(k6/Q)
LET W=IP(k4/Q)
LET kd=IP(U/7+7)
LET kp=IP(SQR(k6))

DIM D(0 TO U-W)

LET COUNT=0
FOR r=1 TO 8
   LET rr=B(r)
   MAT D = ZER
   LET MD=0
   LET t=4

   DO
      LET cca=0
      LET x=G(t)
      IF x^2>k6 THEN EXIT DO
      LET G1=INT(W/x)
      IF MOD(x-rr,Q)=0 THEN
         LET y=(x-rr)/Q
         GOTO 800
      ELSE

         FOR i=1 TO 7
            LET rv=B(i)
            IF MOD(x*rv+rr,Q)=0 THEN
               LET y=-(x*rv+rr)/Q
               GOTO 800
               EXIT FOR
            END IF
         NEXT i
      END IF

800       IF x*G1+y < W THEN
             DO
                LET G1=G1+1
                IF x*G1+y => W THEN EXIT DO
             LOOP
          END IF

          FOR f=G1 TO kd
             IF x*f+y < W THEN  GOTO 900
             IF x*f+y>U THEN GOTO 1000
             LET D(x*f+y-W)=1
900       NEXT f
1000       LET t=t+1

        LOOP

        FOR n=0 TO U-W
           LET ST=n+W
           IF D(n)=0 THEN
              IF Q*ST+rr>k4 AND Q*ST+rr<k6 THEN LET COUNT=COUNT+1
           END IF
        NEXT n
        LET L=L+6
     NEXT r
     LET cprime=COUNT
  END FUNCTION

  EXTERNAL SUB prime(k)
     OPTION ARITHMETIC NATIVE
     DECLARE NUMERIC G(1536277) !素数
     !エラトステネスの篩
     LET Fu=5633
     LET Fm=739
     DIM P(Fu)
     DIM A(Fm)
     MAT P=ZER
     MAT A=ZER
     LET A(1)=2
     LET H1=1
     FOR I=3 TO SQR(Fu) STEP 2
        IF P(I)=0 THEN
           FOR J=I*I TO Fu STEP I
              LET P(J)=1
           NEXT J
        END IF
     NEXT I
     FOR I=3 TO Fu STEP 2
        IF P(I)=0 THEN
           LET H1=H1+1
           LET A(H1)=I
        END IF
     NEXT I

     LET Q=6
     LET k7=k          !篩の計算範囲
     LET k5=IP(k7/Q)+1
     DIM Au(k5),Av(k5)

     MAT Au = ZER     !(6*n-1)
     MAT Av = ZER     !(6*n+1)

     FOR n=3 TO Fm
        LET Pu=A(n)
        IF Pu^2>=k THEN EXIT FOR
        IF MOD(Pu+1,Q)=0 THEN !(6*n-1)
           LET ru=(Pu+1)/Q
           FOR i=1 TO k5
              IF Pu*i+ru>k5 THEN EXIT FOR
              LET Au(Pu*i+ru)=1
           NEXT i
        END IF

        IF MOD(Pu-1,Q)=0 THEN
           LET ru=(Pu-1)/Q
           FOR i=1 TO k5
              IF Pu*i-ru>k5 THEN EXIT FOR
              LET Au(Pu*i-ru)=1
           NEXT i
        END IF

        IF MOD(Pu+1,Q)=0 THEN !(6*n+1)
           LET ru=(Pu+1)/Q
           FOR i=1 TO k5
              IF Pu*i-ru>k5 THEN EXIT FOR
              LET Av(Pu*i-ru)=1
           NEXT i
        END IF

        IF MOD(Pu-1,Q)=0 THEN
           LET ru=(Pu-1)/Q
           FOR i=1 TO k5
              IF Pu*i+ru>k5 THEN EXIT FOR
              LET Av(Pu*i+ru)=1
           NEXT i
        END IF
     NEXT n

     LET G(1)=2
     LET G(2)=3
     LET cz=2
     FOR n=1 TO k5
        IF 6*n-1>k7 THEN GOTO 100
        IF Au(n)=0 THEN
           LET cz=cz+1
           LET G(cz)=6*n-1
        END IF
100    IF 6*n+1>k7 THEN EXIT FOR
       IF Av(n)=0  THEN
          LET cz=cz+1
          LET G(cz)=6*n+1
       END IF
    NEXT n
    PUT TO sha FROM G
    SIGNAL Ok5
END SUB


初期設定は LET E=1E10 !    (1E10)455052511

計算結果

2  11078937  5317482    .156秒
3  16252325  5173388    .141秒
4  21336326  5084001    .152秒
5  26355867  5019541    .142秒
6  31324703  4968836    .168秒
7  36252931  4928228    .145秒
8  41146179  4893248    .163秒
9  46009215  4863036    .164秒
10  50847534  4838319    .165秒

途中省略

90  411523195  4363605    .167秒
91  415885628  4362433    .162秒
92  420243162  4357534    .170秒
93  424603409  4360247    .170秒
94  428958595  4355186    .162秒
95  433311792  4353197    .163秒
96  437663672  4351880    .179秒
97  442014876  4351204    .168秒
98  446362736  4347860    .170秒
99  450708777  4346041    .160秒
100  455052511  4343734    .175秒
   16.48秒     Lazarus fpc-3.2.2-win32

----------------------------------------

2  11078937  5317482    .082秒
3  16252325  5173388    .050秒
4  21336326  5084001    .067秒
5  26355867  5019541    .061秒
6  31324703  4968836    .070秒
7  36252931  4928228    .072秒
8  41146179  4893248    .083秒
9  46009215  4863036    .052秒
10  50847534  4838319    .071秒

途中省略

90  411523195  4363605    .079秒
91  415885628  4362433    .062秒
92  420243162  4357534    .088秒
93  424603409  4360247    .054秒
94  428958595  4355186    .082秒
95  433311792  4353197    .085秒
96  437663672  4351880    .076秒
97  442014876  4351204    .055秒
98  446362736  4347860    .086秒
99  450708777  4346041    .086秒
100  455052511  4343734    .065秒
    7.25秒     Lazarus fpc-3.2.2-win64

動作環境

Intel Core i9 -11900K 自作パソコン

Windows Version
Microsoft Windows 10 (10.0) Professional 64-bit

Paract BASIC Ver. 2.1.2.4 Rev.2   (2021.06.19)

Lazarus fpc-3.2.2-win64
lazarus-2.2.0RC1-fpc-3.2.2-win64.exe               2021-07-08 199.7 MB
lazarus-2.2.0RC1-fpc-3.2.2-cross-i386-win32-win64.exe 2021-07-08  55.1 MB
 

レンタル掲示板
/73