***Introducció de valors per realitzar el primer intent de sinus
#define Pi 3.141592653589793238463
***Interval de valors i els seus sinus corresponents
*#define cX {0.0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,1.1,1.2,1.3,1.4,1.5,1.6,1.7}
**#define cY {0.0,0.099833416646828115231, 0.19866933079506121545,;
*            0.29552020666133957510,0.38941834230865049166,;
*            0.47942553860420300027,0.56464247339503533720,;
*            0.64421768723769105367,0.71735509089952276162,;
*            0.78332690962748338846,0.84147098480789650665,;
*            0.89120736006143533995,0.93203908596722634967,;
*            0.96355818541719296470,0.98544972998846018065,;
*            0.99749498660405443094,0.99957360304150516434,;
*            0.9916648105}
#define nCorba 3
*Definició de variables que seran constants durant tot el programa
STATIC Pt,Vix,Viy,xi,yi,dt,m,ax,ay,Fx,Fy,Imprimir,Sortir,Zoom,;
       Hescala,Emec
*******************************************************************************
Main()
*******************************************************************************
PROCEDURE Main
Zoom:=100
r:=SetGMode(6)
*Fixació de valors recomanables per un exemple determinat
Hescala:=0.1  && 0.10 unitats per pixel
Sortir:="N"
*Imprimir:="N"
Viy:=-15
Vix:=15
xi:=25
yi:=25
dt:=0.002
m:= 1
Fx:="0                                           "+space(150)
Fy:="0                                           "+space(150)
*Pt:=space(150)
Pt:="-51000/(x^2+y^2)^0.5                                 "
*Ordre d'operacions del programa.
*Pantalla de presentació
Presentacio()

*Ordre per sortir
  DO WHILE Sortir=="N"
    Valors()
*L'ordinador passa a mode gràfic
    r=SetGMode(18)
*L'ordinador passa a alta resolució
    r=SetHiRes(0)
*Ordre per decidir que si hi ha alguna dada escrita en les forces que faci
*els calculs de fFemfer
      IF LEN(ALLTRIM(Fx))>0 .AND. LEN(ALLTRIM(Fy))>0
*Realització del progama depenent de la força
        r=fFemFer()
      ELSE
*Realització del programa depenent de l'energia potencial
        r=vFemFer()
      ENDIF
*Ordre per tornar l'ordinador a mode text
    r=SetText()
*Rutines de sortida:
*Realització de l'ordre vull sortir
    VullSortir()
  ENDDO
*Pantalla de sortida
pntsortida()
*Neteja de la pantalla
CLEAR SCREEN
RETURN
*******************************************************************************
*Realització de la pantalla de portada
FUNCTION presentacio()
*Neteja la pantalla
CLEAR SCREEN
r=SetGMode(18)
r=SetHiRes(0)
r=PicRead(100,50,0,"C:camps2.pcx")
*Ordres per música d'entrada (correspon a les freqüències de cada nota)
TONE(440,6)
TONE(440*2^(4/12),6)
TONE(440*2^(2/12),6)
TONE(440,6)
TONE(440*2^(-1/12),6)
TONE(440*2^(-3/12),6)
TONE(440*2^(-1/12),6)
TONE(440,16)
r=SetText()
RETURN NIL
******************************************************************************
*Taula per introduir els valors del programa
FUNCTION Valors()
CLEAR SCREEN
*Ordres per fer aparèixer els enunciats per introduir els valors corresponents
@ 0,2 SAY "INTRODUIU ELS VALORS DESITJATS:" COLOR "R"
@ 1,0 TO 23,75 DOUBLE COLOR "R+/B+"
@ 2,2 SAY "Posició inicial en l'eix de les 'x': " GET xi PICTURE "999999.999999"
@ 4,2 SAY "Posició inicial en l'eix de les 'y': " GET yi PICTURE "999999.999999"
@ 6,2 SAY "Velocitat en l'eix de les 'x': " GET Vix PICTURE "999999.999999"
@ 8,2 SAY "Velocitat en l'eix de les 'y': " GET Viy PICTURE "999999.999999"
@ 10,2 SAY "Massa de la partícula: " GET m PICTURE "999999.999999"
@ 12,2 SAY "Força a l'eix de les 'x': " GET Fx PICTURE "@S44"
@ 14,2 SAY "Força a l'eix de les 'y': " GET Fy PICTURE "@S44"
@ 16,2 SAY "Fórmula de l'Energia Potencial:" GET Pt PICTURE "@S39"
@ 18,2 SAY "Diferencial de temps: " GET dt PICTURE "999999.9999999"
@ 20,2 SAY "Separació entre unitats: " GET Zoom PICTURE "9999"
@ 22,2 SAY "Escala de la pantalla (unitat per píxel): " GET Hescala PICTURE "9999.9999"
READ
RETURN NIL
******************************************************************************
*Funció per captar la fórmula de la força en l'eix "x"
FUNCTION cFx(x,y,vx,vy,t)
LOCAL fFormula,fCodi
fFormula:="{|x,y,vx,vy,t| "+ALLTRIM(Fx)+"}"
fCodi:=&(fFormula)
RETURN EVAL(fCodi,x,y,vx,vy,t)
*******************************************************************************
*Funció per captar la fórmula de la força en l'eix "y"
FUNCTION cFy(x,y,vx,vy,t)
LOCAL fFormula,fCodi
fFormula:="{|x,y,vx,vy,t| "+ALLTRIM(Fy)+"}"
fCodi:=&(fFormula)
RETURN EVAL(fCodi,x,y,vx,vy,t)
******************************************************************************
*Realització dels càlculs de la posició de la partícula a partir de les forces
*(aquí els procesos inicials són iguals que a l'apartat anterior)
FUNCTION fFemFer()
LOCAL n,r,t,ax,ay,bx,by,Lambda,Emec,Ecin,dl,cCordX,cCordY
dl:=0.0000000000001
Eixos()
n:=0
xiVell:=xi
yiVell:=yi
  FOR t:=0 TO 1500 STEP dt
      IF n%10==0
        cCordX:=xiVell/Hescala+650
        cCordY:=yiVell/Hescala+500
        r:=DrawCircle(cCordX,cCordY,3,0,360,8,0,0)
        r:=SetPix(cCordX,cCordY,0,15)
      ENDIF
*Càlcul de les acceleracions a partir de la força i de la massa
    ax:=cFx(xi,yi,Vix,Viy,dt)/m
    ay:=cFy(xi,yi,Vix,Viy,dt)/m
*Fórmules de cinemàtica per calcular la posició i la velocitat
    xi:=xi+Vix*dt+(ax*dt^2)/2
    yi:=yi+Viy*dt+(ay*dt^2)/2
    Vix:=Vix+ax*dt
    Viy:=Viy+ay*dt
      IF n%10==0
*Ordre per fer la representació de la velocitat centrada a la pantalla.
        r:=SetPix(Vix+650,Viy+500,0,14)
        cCordX:=xi/Hescala+650
        cCordY:=yi/Hescala+500
*Ordre per tal que quan la partícula arribi als marges de la pantalla no la dibuixi
          IF cCordX<1300 .AND. cCordX>-1 .AND. cCordY<1000 .AND. cCordY>-1
*Ordre per dibuixar la partícula
            r:=DrawCircle(cCordX,cCordY,3,0,360,8,0,12)
*Fi de les ordres.
          ENDIF
*Reconeixement de les posicions per fer els nous càlculs
        xiVell:=xi
        yiVell:=yi
        INKEY(0.1)
          IF LASTKEY()==27
            RETURN
          ENDIF
      ENDIF
    n:=n+1
  NEXT t
INKEY(0)
RETURN
*******************************************************************************
*Funció per captar l'energia potencial
FUNCTION Potencial(x,y)
LOCAL fFormula,fCodi
fFormula:="{|x,y| "+ALLTRIM(Pt)+"}"
fCodi:=&(fFormula)
RETURN EVAL(fCodi,x,y)
******************************************************************************
*Subapartat de càlcul de la posició i la velocitat a partir de l'E.potencial
FUNCTION vFemFer()
*Definició de variables que actuaran en aquest apartat
LOCAL n,r,t,ax,ay,bx,by,Lambda,Emec,Ecin,dl,cCordX,cCordY
*Diferencial de temps imposat per fer els càlculs de la derivada
dl:=0.00000000001
*Apartat de càlcul
Eixos()
*Fórmula de l'energia mecànica a partir de la potencial i la cinètica, amb
*les velocitats en l'eix de les x i de les y
Emec:=Potencial(xi,yi)+(m*(Vix^2+Viy^2))/2
n:=0
*Definició de les variables de la posició
xiVell:=xi
yiVell:=yi
*Límit de temps
  FOR t:=0 TO 1500 STEP dt
*Ordre per tal que actui cada vegada que t sigui multiple de 10
      IF n%10==0
*Definició de la posició de la partícula centrant-la a la pantalla (la nostra
*pantalla té 1000x1500 pixels) d'aquí que el centre és:(650,500)
        cCordX:=xiVell/Hescala+650
        cCordY:=yiVell/Hescala+500
*Ordre per borrar el dibuix de la partícula abans de ser dibuixada per primera
*vegada
        r:=DrawCircle(cCordX,cCordY,3,0,360,8,0,0)
*Ordre per fer el rastre de la partícula (en punts blancs)
        r:=SetPix(cCordX,cCordY,0,15)
      ENDIF
*Definició de l'acceleració en l'eix x i a l'eix y a partir de la deriva de
*l'energia potencial
    ax:=-(Potencial(xi+dl,yi)-Potencial(xi,yi))/(m*dl)
    ay:=-(Potencial(xi,yi+dl)-Potencial(xi,yi))/(m*dl)
*    xi:=xi+Vix*dt+(ax*dt$\hat{\phantom{n}}$2)/2
*    yi:=yi+Viy*dt+(ay*dt$\hat{\phantom{n}}$2)/2
*Fórmules del moviment uniforme
    xi:=xi+Vix*dt
    yi:=yi+Viy*dt
*Definició de l'energia cinètica a partir de la resta de l'energia mecànica
*i de la potencial
    Ecin:=Emec-Potencial(xi,yi)
*Càlcul de les noves acceleracions per obtenir les noves velocitat
    ax:=-(Potencial(xi+dl,yi)-Potencial(xi,yi))/(m*dl)
    ay:=-(Potencial(xi,yi+dl)-Potencial(xi,yi))/(m*dl)
*Fórmules de la velocitat en Moviment uniformement accelerat dels nous punts
    bx:=Vix+ax*dt
    by:=Viy+ay*dt
*Coeficient de restitució de la pèrdua d' energia
    Lambda:=(((Emec-Potencial(xi,yi))*2)/(m*((Vix+(ax*dt))$\hat{\phantom{n}}$2+(Viy+(ay*dt))$\hat{\phantom{n}}$2)))$\hat{\phantom{n}}$(1/2)
*Aplicació del valor Lambda sobre les velocitats
    Vix:=Lambda*bx
    Viy:=Lambda*by
*Representació gràfica de la velocitat de la partícula en cada instant
      IF n%10==0
        r:=SetPix(Vix+650,Viy+500,0,14)
*Posició de la velocitat en l'eix de coordenades
        cCordX:=xi/Hescala+650
        cCordY:=yi/Hescala+500
*Instruccions per restringir la representació gràfica dins la pantalla
          IF cCordX<1300 .AND. cCordX>-1 .AND. cCordY<1000 .AND. cCordY>-1
*Representació de la partícula
            r:=DrawCircle(cCordX,cCordY,3,0,360,8,0,4)
          ENDIF
*Reconeixement de les posicions per fer de nou els càlculs
        xiVell:=xi
        yiVell:=yi
        INKEY(0.1)
          IF LASTKEY()==27
            RETURN
          ENDIF
      ENDIF
    n:=n+1
  NEXT t
INKEY(0)
RETURN
*******************************************************************************
*Funció per tal d'establir els eixos de coordenades
FUNCTION Eixos()
*En aquest apartat es dibuixa la part negativa dels eixos y
  FOR e:=500 TO 0 STEP -Zoom
    SayString(660,e,1,0,13,ALLTRIM(STR((e-500)*Hescala,5,0)))
    DrawLine(0,e,1500,e,0,1,11)
  NEXT e
*En aquest apartat es dibuixa la part positiva dels eixos y
  FOR e:=500 TO 1000 STEP Zoom
    SayString(660,e,1,0,13,ALLTRIM(STR((e-500)*Hescala,5,0)))
    DrawLine(0,e,1500,e,0,1,11)
  NEXT e
*En aquest apartat es dibuixa la part negativa dels eixos x
  FOR e:=650 TO 0 STEP -Zoom
    SayString(e+10,500,1,0,13,ALLTRIM(STR((e-650)*Hescala,5,0)))
    DrawLine(e,0,e,1000,0,1,11)
  NEXT e
*En aquest apartat es dibuixa la part positiva dels eixos x
  FOR e:=650 TO 1500 STEP Zoom
    SayString(e+10,500,1,0,13,ALLTRIM(STR((e-650)*Hescala,5,0)))
    DrawLine(e,0,e,1000,0,1,11)
   NEXT e
*Dibuix dels eixos principals
DrawLine(650,1000,650,0,0,0,9)
DrawLine(0,500,1500,500,0,0,9)
RETURN NIL
******************************************************************************
*Funció per preguntar si es vol continuar o si contrariament es vol sortir
FUNCTION VullSortir()
CLEAR SCREEN
@ 10,5 SAY "VOLEU SORTIR? (S/N)" GET Sortir
READ
Sortir:=UPPER(Sortir)
RETURN NIL
******************************************************************************
*Primer intent fracassat per imprimir les pantalles de representació, cal destacar
*que posteriorment vam utilitzar el Capturador d'imatges del WP.5.1
*FUNCTION VullImprimir()
*CLEAR SCREEN
*@ 5,5 Say "VOLEU IMPRIMIR? (S/N)" GET Imprimir
*READ
*Imprimir:=UPPER(Imprimir)
*  IF Imprimir=="S"
*    r=PrintPCL(2,100,200,600)
*  ENDIF
*RETURN NIl
******************************************************************************
*Realització de la pantalla de sortida
FUNCTION pntsortida()
*Neteja la pantalla
CLEAR SCREEN
r=SetGMode(18)
r=SetHiRes(0)
r=PicRead(100,110,0,"C:camps1.pcx")
*Ordres per música de sortida (correspon a les freqüències de cada nota)
TONE(440,6)
TONE(440*2^(4/12),6)
TONE(440*2^(2/12),6)
TONE(440,6)
TONE(440*2^(-1/12),6)
TONE(440*2^(-3/12),6)
TONE(440*2^(-1/12),6)
TONE(440,16)
r=SetText()
RETURN NIL
******************************************************************************
*Funció per calcular el sinus
FUNCTION sin(aPidpi)
*Definició de variables que intervindran en aquest apartat
LOCAL sSigne,valSinus,onSoc,nNum,nDen,nI,nJ
sSigne:=1
*Si el valor del sinus és més petit que 0
  IF aPidpi<0
*llavors, el canvia de signe
    aPidpi:=-aPidpi
    sSigne:=-1
  ENDIF
*Reducció del valor introduït entre 0 i 2*pi, que és la part on definirem
*la funció sinus amb el polinomi de grau onze
aPidpi:=aPidpi%(2*Pi)
*Si el valor introduït és, encara, més gran que pi li resta pi
  IF aPidpi>Pi
    aPidpi:=aPidpi-Pi
*I si queda més petit que 0 el canvia de signe
    sSigne:=-sSigne
  ENDIF
*Si el sinus és més gran que pi/2
  IF aPidpi>Pi/2
*llavors, que el valor del sinus sigui pi menys el valor introduït
    aPidpi:=Pi-aPidpi
  ENDIF
*Fórmula del polinomi d'onzè grau per calcular el sinus
RETURN sSigne*(-(aPidpi^11)/(2*3*4*5*6*7*8*9*10*11);
              +(aPidpi^9)/(2*3*4*5*6*7*8*9);
              -(aPidpi^7)/(2*3*4*5*6*7);
              +(aPidpi^5)/(2*3*4*5);
              -(aPidpi^3)/(2*3);
              +aPidpi)
*******************************************************************************
*Funció per calcular el cosinus (a partir del sinus)
FUNCTION cos(aPidpi)
RETURN sin(aPidpi+(Pi/2))
*******************************************************************************
*Funció per calcular la tangent (a partir del sinus i del cosinus)
FUNCTION tan(aPidpi)
RETURN sin(aPidpi)/cos(aPidpi)
*Casella per comprobar l'efectivitat del sinus, cosinus i tangent. Actualment
*esborrades.
/*
LOCAL Hangle
Hangle:=0
CLEAR SCREEN
@ 2,2 GET Hangle PICTURE "9999999999.99999999999999"
READ
@ 4,2 SAY "sin ="+STR(sin(Hangle),17,14)
@ 5,2 SAY "cos ="+STR(cos(Hangle),17,14)
@ 6,2 SAY "tan ="+STR(tan(Hangle),17,14)
*/
******************************************************************************