10 REM OSCILLATORE SMORZATO INTEGRATO CON RUNGE KUTTA 4 - FULL
20 REM DATA: PASSI, DURATA S, FREQUENZA HZ, SMORZAMENTO
30 DATA 400,2,1,0.1
40 LET NSTEPS=0:LET TFINAL=0:LET FREQ=0:LET DAMP=0
50 READ NSTEPS,TFINAL,FREQ,DAMP
60 IF NSTEPS<20 OR NSTEPS>2000 THEN GOTO 510
70 IF NSTEPS<>INT(NSTEPS) THEN GOTO 510
80 IF TFINAL<=0 OR FREQ<=0 THEN GOTO 510
90 IF DAMP<0 OR DAMP>=1 THEN GOTO 510
100 LET OMEGA=2*3.141592653589793*FREQ
110 LET DT=TFINAL/NSTEPS
120 IF DT*OMEGA>0.1 THEN GOTO 510
130 LET POS=1:LET VEL=0:LET TIME=0
140 LET XP=0:LET VP=0:LET DX=0:LET DV=0
150 LET KONE=0:LET LTWO=0:LET LONE=0:LET KTWO=0
160 LET KTHREE=0:LET LTHREE=0:LET KFOUR=0:LET LFOUR=0
170 LET EXACT=0:LET ERRMAX=0:LET ENERGY=0:LET OLDEN=0
180 LET ALPHA=DAMP*OMEGA:LET BETA=OMEGA*SQR(1-DAMP*DAMP)
190 LET STRIDE=INT(NSTEPS/10):LET COUNT=0
200 OLDEN=OMEGA*OMEGA/2
210 FOR I=1 TO NSTEPS
220 XP=POS:VP=VEL:GOSUB 450
230 KONE=DX:LONE=DV
240 XP=POS+DT*KONE/2:VP=VEL+DT*LONE/2:GOSUB 450
250 KTWO=DX:LTWO=DV
260 XP=POS+DT*KTWO/2:VP=VEL+DT*LTWO/2:GOSUB 450
270 KTHREE=DX:LTHREE=DV
280 XP=POS+DT*KTHREE:VP=VEL+DT*LTHREE:GOSUB 450
290 KFOUR=DX:LFOUR=DV
300 POS=POS+DT*(KONE+2*KTWO+2*KTHREE+KFOUR)/6
310 VEL=VEL+DT*(LONE+2*LTWO+2*LTHREE+LFOUR)/6
320 TIME=I*DT
330 EXACT=EXP(-ALPHA*TIME)*(COS(BETA*TIME)+ALPHA*SIN(BETA*TIME)/BETA)
340 IF ABS(POS-EXACT)>ERRMAX THEN ERRMAX=ABS(POS-EXACT)
350 ENERGY=(VEL*VEL+OMEGA*OMEGA*POS*POS)/2
360 IF ENERGY>OLDEN+0.000001 THEN GOTO 520
370 OLDEN=ENERGY:COUNT=COUNT+1
380 IF I-STRIDE*INT(I/STRIDE)=0 THEN GOSUB 490
390 NEXT I
400 PRINT "FINAL;",TIME,";",POS,";",VEL,";",EXACT,CRLF
410 PRINT "MAXERR;",ERRMAX,";COUNT;",COUNT,CRLF
420 IF ERRMAX>0.0001 OR COUNT<>NSTEPS THEN GOTO 520
430 PRINT "RK4 TEST OK",CRLF
440 END
450 REM DERIVATE DELLO STATO: DX=VELOCITA, DV=ACCELERAZIONE
460 DX=VP
470 DV=-2*DAMP*OMEGA*VP-OMEGA*OMEGA*XP
480 RETURN
490 PRINT "STATE;",TIME,";",POS,";",VEL,";",ENERGY,CRLF
500 RETURN
510 PRINT "ERRORE PARAMETRI RK4",CRLF:END
520 PRINT "ERRORE AUTOVERIFICA RK4",CRLF:END