VPython

[Ivo_Timóteo]

Tal como prometido, vou apresentar o VPython à comunidade Quarkiana.

Vou assumir que já sabem um pouco de Python e passar à frente as partes específicas da linguagem (que acabam por ser muito poucas )

Em primeiro lugar têm de ter o VPython instalado no vosso computador.

Windows Users:
Podem ir ao site do VPython e fazer download da última package estável, mesmo que ainda n tenham o Python instalado. (Ela inclui o
IDLE, o modulo visual e outras coisas... )

Linux Users:
A partida já têm o Python Para terem o módulo visual, consultem a página

Vamos então começar por uma simples simulação de uma bola em queda livre que ressalta ao tocar no chão e tem um coeficiente de restituição igual a 1.


Em primeiro lugar, importamos tudo do módulo visual:

Código: Seleccionar Todos

from visual import *

Depois definimos umas variáveis que nos vão ser úteis

Código: Seleccionar Todos

g = 9.8
dt = 0.01

Agora criamos dois objectos, um vai representar a bola que é, naturalmente, uma esfera. O chão vai ser uma caixa.
Os argumentos parecem-me bastante óbvios mas nada como dar aqui um saltinho para ficar a saber tudo sobre como
se inicializam e como interagem os objectos e janelas no VPython.

Código: Seleccionar Todos

bola = sphere(pos=(0., 20., 0.), radius=2.0, color=color.red)
terra = box(pos=(0., -0.5, 0.), size=(20., 0.5, 20.), color=color.blue)

Neste momento estou a definir a variável vel associada ao objecto bola. Naturalmente vai representar a velocidade da bola.
O tipo de dados vector está definido em visual e dá bastante jeito, vejam na reference. Para além de ter várias operações
sobre vectores definidas, pode-se aceder às diferentes coordenadas com o nome x, y e z. Ex. bola.vel.y

Código: Seleccionar Todos

bola.vel = vector(0.,0., 0.)

Aqui controlo a "janela" que tem como nome default "scene"

Código: Seleccionar Todos

scene.autoscale=0
scene.center=(0,10,0)
scene.range=(25,25,25)

finalmente o ciclo de execução do programa:

Código: Seleccionar Todos

while 1:
    rate(100) # faz com que o programa execute 100 ciclos por segundo
                  # com dt= 0.01, temos "tempo real"
    bola.pos = bola.pos + bola.vel*dt  #esta ja todos conhecemos :)
    if bola.pos.y < bola.radius:   #condicao de colisao
        bola.vel.y = -bola.vel.y
    else:
        bola.vel.y = bola.vel.y -g*dt

Finalmente, se vos parecer demasiado lento para tempo real, lembrem-se que é uma esfera de 4 metros de diâmetro a cair de 20 metros de altura

Comentários, correcções e perguntas são benvindas!

[Ivo_Timóteo]

Agora uma pequena evolução! Vamos introduzir multiplas janelas e comandos

Julgo que o código não deve levantar muitas dúvidas. Não se esqueçam de dar uma vista de olhos à reference do VPython.
Qualquer dúvida estão à vontade

Eu pessoalmente gosto bastante de ver exemplos

Código: Seleccionar Todos

from visual import *
from visual.controls import *

#constantes
g = 9.8
dt = 0.01

#lado da janela
lado = 350

#janela principal
display(x=0, y=0, width=lado, height=lado, range=(25,25,25), center=(0,10,0))
#janela dos controlos
c = controls(x=lado, y=0, width=lado, height=lado, range=60,
             title="coef. restituicao");

#inicializacao dos objectos
bola = sphere(pos=(0., 20., 0.), radius=2.0, color=color.red)
terra = box(pos=(0., -0.5, 0.), size=(20, 0.5, 20), color=color.blue)

#algumas propriedades da bola
bola.vel = vector(0.,-0.1, 0.)
bola.cr = 1

#evita o autoscale
scene.autoscale=0

#funcao que altera o coeficiente de restituicao
def setC_r(x):
    bola.cr=x.value

#reset da posicao da bola
def setPos():
    bola.pos.y = 20.
    bola.vel.y = 0.

#inicializacao do slider
s = slider(pos=(-15, 40), width=10, length=80, axis=(0,-1,0), min=0, max=1,
           action=lambda: setC_r(s))

#inicializacao do reset
b = button(pos=(25, 0), witdth=15, length=20, text="RESET", action=lambda: setPos())

#escolha do valor inicial do slider
s.value=1


while 1:
    #controla o rate do ciclo (deve ser usado com atencao ao dt
    rate(100)
    #procura interaccoes nos controlos
    c.interact()
    #actualiza a posicao
    bola.pos = bola.pos + bola.vel*dt

    #teste de colisao
    if bola.pos.y < bola.radius:
        bola.pos.y=bola.radius
        bola.vel.y = -bola.vel.y*bola.cr
    else:
        bola.vel.y = bola.vel.y -g*dt

[RicardoCampos]

Há VPython para 2.6 no Windows?

[Ivo_Timóteo]

Não sei, vê se a última versão do site é compatível

[jap]

Creio que a última versão disponível para Windows do VPython só funciona (por enquanto) com o Python 2.5. Na realidade, não experimentei ainda usá-lo com o Python 2.6, é preciso testar.

[Ivo_Timóteo]

Nova versão do vpython Visual 5 aqui.

[Miro]

já não pegava em python desde a fisica@UC no verão, decidi pegar ontem e exprimentar o Vpython, isto é bem porreiro!

fiz um programita simples, deve estar um bocado "estouvado" mas cumpre o objectivo

simula um "xuto" numa bola (um bloco vermelho a ir contra uma bola conta como xuto?)

inserimos a altura a que se dá o xuto, a velocidade imprimida pelo xuto, o angulo, aquilo simula a 3D e dá alguns dados, como tempo de voo, distancia etc.

agora tenho e incluir ressaltos, atritos do ar e por ai fora

deixo aqui para se alguem quiser ver

http://www.megaupload.com/?d=04S3UA9B

um screenshot:



fica aqui o codigo (aumentei o tamanho do chao, no programa e na pic está minusculo)

Código: Seleccionar Todos

from visual import *

# tamanho do ecra
scene.width = 800
scene.height = 600

# parte do espaço 3D k mostra

#scene.autoscale = 0
scene.range = (300,200,100)
scene.center = (50,50,0)

valuesset = true

mylabel = label(pos=(100,100,0),text="velocity x: ",height=10)
mylabel.text = "velocity x: 0" + "\nvelocity y: 0 "

y=2
n = 0

#objectos
ball = sphere(pos=(0,4,0),radius=4, color=color.green,material=materials.earth) #uma bola
ground = box(pos=(50,-1,0),size=(1000,0.5,500),material=materials.wood) #o chao
canhao = box(pos=(-70,2,0),size=(10,5,7), color=color.red,material=materials.wood) #o canhao simbolico

#constantes e coiso
gravity = 9.8  # m/s**2

#user settings

if valuesset == true:
    velocity=input('Introduz um valor para a velocidade imprimida pelo "xuto" (m/s) :')
    angle=input('Introduz um valor para o angulo do "xuto" (graus) :')
    inicialalt=input('A que Altura do chão se dá o "xuto"? (m) : ')
    altlaicnh=inicialalt
    inicialalt= inicialalt + 4
    w=input('queres que a bola deixe uma trajectoria? y ou n : ')
    valuesset = false
   
#inicialalt=4
#angle = 50
#velocity = 20
satan =angle #pa no fim nao ter de converter pa graus outra bez
angle = angle * (pi/180) # convertido pa radianos

#ajustar altura inicial
canhao.pos.y = inicialalt
ball.pos.y = inicialalt

VelocityY = velocity * sin(angle)
VelocityX = velocity * cos(angle)

seconds = 0
dt = .01 #incremento

if not valuesset:

    #loop do canhao
    disparou = False
    while not disparou:
        rate(100)
        seconds += dt
        canhaoX = -70 + velocity * seconds
        canhao.pos = vector(canhaoX,inicialalt,0) #pa mudar canhao de sitio no coiso 3D

        if canhao.pos.x >= 0:
             disparou = True
             seconds = 0
             
        if w >= 1:
            ball.trail = curve(color=color.cyan)#rasto
       

    #loop da bola

    finished = False
    while not finished:
        rate(100) # nao correr mais de 100 vezes per sec
        seconds += dt

        # y(t) = y0 + v0*t + .5 * a * t**2
        # v(t) = v0 + a * t
        ballY = inicialalt + VelocityY * seconds - .5 * gravity * seconds**2
        ballX = VelocityX * seconds
        v = VelocityY - gravity * seconds

        if VelocityX <=0.005:
            if VelocityX >= -0.005:
                VelocityX = 0

        ball.pos = vector(ballX,ballY,0) #pa mudar a bola de sitio no coiso 3D

        if w >= 1 :
            ball.trail.append(pos=ball.pos) #adicionar rast0

        if v >= 0: #se a velocidade for positiva
            alturamax = ballY #altura maxima é igual á altura a que a bola tá
            #ao deixar de ser positiva deixa de mudar e sabemos a altura atingida
       
        mylabel.text = "velocity x: " + str(VelocityX) + "\nvelocity y: " + str(v)


        velocidadesolo = - v
   
        if ballY - 4 <= 0: #o menos 4 é o raio da bola, senao fazia as contas para centro d bola
            finished = True #acaba o loop
            mylabel.text = "velocity x: 0" + "\nvelocity y: 0 "
            if ballX <= 0.005: #obtinha valores estranhos sem isto tipo 3.12580771813e-014
                if ballX >= -0.005:
                    ballX = 0
            # agora diz estas tretas todas

            data = label(pos=(-70,80,0),text="data",height=10)
            data.text = "velocidade inicial: " + str(velocity) + "\nangulo lancamento: " + str(satan) + "\ntempo de voo: " + str(seconds) + "\ndistancia percorrida segundo X : " + str(ballX) + "\naltura maxima atingida : " + str(alturamax) + "\nvelocidade com que atingiu o solo segundo y : " +str(velocidadesolo)       
       


[Miro]

apos adicionar 3 linhas o meu programa deixou de dar direito, os objectos a certa altura "encolhem" todos ate so ver um ecran preto

ja apaguei essas linhas extra de codigo e nada....ja reiniciei o python varias vezes e tudo

ele estava a dar perfeitamente!

(desculpem o double post)

[Bruno Oliveira]

Se calhar podias pôr aqui o código e algum VPython expert pudesse ajudar-te no debug dele.

[Miro]

bem, já resolvi o problema

a velocidade é dado pela equação v = v0 + at

se eu meter o v0 acontece isso mas se o tirar não.....

vou masé fazer o programa de raiz, estava a reciclar coisas de outro lado deve estar ai o problema....

fica ai o codigo anyways (nao liguem aos nomes parvos )

Código: Seleccionar Todos

from visual import *

# propriedades da janela
scene.width = 800
scene.height = 600

scene.autoscale = 0 # para o vpython nao zoomar sozinho á toa
scene.range = (100,100,100) #quanto mostra na janela
scene.center = (50,40,0) # centra camera neste ponto

# objecto e valores
bixo = cylinder(pos=(0,0,0),axis=(0,2,0), radius=2, color=color.green)
bixo.velocity = vector(0,0,0)
bixo.acceleration = vector(0,0,0)
bixo.mass = 10 # kg
bixo.pos = vector (0,0,0)

# criar um xão
ground = box(pos=(50,-1,1),size=(500,2,500))


# vectores de aceleração
accelerationGravity = vector(0,-9.8,0)  # m/s**2
accelerationNormal = vector(0,0,0) # m/s**2,
accelerationKineticFriction = vector(0,0,0) # m/s**2, ate xegar ao xão não tem disto

# condições de lançamento
velocityThrownMagnitude=input('Introduz um valor para a velocidade imprimida pelo "xuto" (m/s) :')
angleThrown=input('Introduz um valor para o angulo do "xuto" (graus) :')
inicialalt=input('A que Altura do chão se dá o "xuto"? (m) : ')
y = 5
n = 0
w=input('queres que o corpo deixe uma trajectoria? y ou n : ')
satan=angleThrown #pa depois indicar antes de tar em rads
angleThrown = angleThrown * (pi/180) # em radianos
bixo.pos= vector (0,inicialalt,0)

#trajectoria?
if w >= 1:
    bixo.trail = curve(color=color.cyan)#rasto

# converter a velocidade de lançamento num vector
velocityThrown = vector(0,0,0)
velocityThrown.y = velocityThrownMagnitude * sin(angleThrown)
velocityThrown.x = velocityThrownMagnitude * cos(angleThrown)
velocityThrown.z = 0

# velocidade do vento :D
print "insira a velocidade do vento segundo cada uma das coordenadas (m/s)"
q = input ('segundo x')
w = input ('segundo y')
e = input ('segundo z')
velocityWind = vector(q,w,e)

# velocidade e aceleraçao do bixo
bixo.velocity = velocityThrown + velocityWind
bixo.acceleration = accelerationGravity + accelerationNormal + accelerationKineticFriction

seconds = 0 #  tempo começa ao zero
dt = .01   

#loop
finished = False
while not finished:
    rate(100) #nao fazer isto mais do que 100 bezes/s
    seconds += dt  # tempo total

   
    #equação da velocidade: v = v0 + at
    # equaçao posiçao: pos = pos0 + v0t + .5 * s * t**2
    #bixo.velocity += velocityThrown + bixo.acceleration * dt #esta e que dá erro
    bixo.velocity += bixo.acceleration * dt
    bixo.pos += bixo.velocity * dt + .5 * bixo.acceleration * dt**2

    #adicionar pontos á curva k forma a trajectoria
    if w >= 1 :
        bixo.trail.append(pos=bixo.pos)
   
           
    #atingiu xao tem atrito
    if bixo.pos.y <= 0:
   

         bixo.velocity.y = 0
         bixo.pos.y = 0
       
        # intensidade da força devido a Fricção = força normal * CoefficientKinetic
         CoefficientOfKineticFriction = .45

        # segunda lei Newton = F=ma
         forceGravity = bixo.mass * accelerationGravity

        # terceira lei de newton
         forceNormal = -forceGravity

        # 2 lei
         accelerationNormal = forceNormal / bixo.mass

        # magnitude da força de fricção
         forceKineticFrictionMagnitude = mag(forceNormal * CoefficientOfKineticFriction)

        # direção da força de fricçao
         forceKineticFrictionDirection = norm(bixo.velocity) * -1

        # multiplica pa ter vector
         forceKineticFriction = forceKineticFrictionMagnitude * forceKineticFrictionDirection

        # arranjar aceleração       
         accelerationKineticFriction = forceKineticFriction / bixo.mass

        # updapte aceleraçao do bixo
         bixo.acceleration = accelerationNormal + accelerationGravity + accelerationKineticFriction
   
    # acabar loop quando para de mexer
    # devido a imprecisões vamos considerar perto de zero
    if mag(bixo.velocity) <= .02:
        finished = True
        bixo.velocity = (0,0,0)

    # mostrar cenas numa janela

data = label(pos=(-70,80,0),text="data",height=10)
data.text ="velocidade inicial: " + str(velocityThrownMagnitude) + "\nangulo lancamento: " + str(satan) + "\ntempo de voo: " + str(seconds) + "\ndistancia percorrida segundo X : " + str(bixo.pos.x)       



     



[jap]

bem, já resolvi o problema

a velocidade é dado pela equação v = v0 + at

se eu meter o v0 acontece isso mas se o tirar não.....
(...)

Há uma muito boa razão para o teu programa não funcionar com a linha onde implementas v = v0 + at e que agora está comentada no código!

Quem é o primeiro a explicar a razão desse bug?

[Tharis]

Li assim o código por alto e tentei perceber por alto o que faz. Supostamente temos um objecto que é atirado com velocidade inicial v0 e que depois tem a aceleração gravítica e o vento ao barulho.

Se bixo.velocity é a velocidade de bixo, então, se ja foi inicializada a variável com o valor de v0, porque é que se havia de estar a cada dt a adicionar-lhe de novo?

Será isto ou não percebi o programa?


Já agora Miro, quando encontrar o link de uma espécie de "manual de bom Python Coding" aconselho-te a ler.

EDIT

Ora aqui está: PEP 8

[jap]

(...)
Se bixo.velocity é a velocidade de bixo, então, se ja foi inicializada a variável com o valor de v0, porque é que se havia de estar a cada dt a adicionar-lhe de novo?

Será isto ou não percebi o programa?
(...)

É isso mesmo!

Por isso é que a linha substituta de código ... está absolutamente correcta, e a supostamente correcta está errada.

Ou seja:

no ciclo apenas implementamos o incremento da velocidade,

v += a*dt

a inicialização de v (atribuindo a v o valor da velocidade inicial) só deve ser feito uma vez e não dentro do ciclo!

Repara que

v + = a*dt significa, na realidade

v_novo = v_antigo + a*dt ...

onde o v_antigo é a velocidade do ciclo anterior!

Deu para perceber, Miro?

[Miro]

oh, estupidez a minha! estava a ignorar completamente o + antes do =

sim já entendi! obrigado pela explicação professor e Tharis

vou ler esse PEP quando tiver tempo, eu sempre ui desorganizado em tudo mas alguma coisa hade adiantar

[Tharis]

Fiz esta coisinha como uma preparação para o hyper-mega projecto que hei-de fazer conjuntamente com people daqui.

O código tá um bocadinho forçado e tá assim à papo-seco as coisas lá.... Utilizei a 1a lei de kepler e pouco mais de interessante...

Código: Seleccionar Todos

#!/usr/bin/env python2.6.2
#
#  Solar.py v0.1 (following the rule "if it doesn't work, name it v0.1")
# 
#
#  Created by Francisco Huhn on 6/8/09
#  Should mention Satelite.py by Joao Casimiro
#  Copyright (c) 2009. All rights reserved.
#

""" This is a script that simulates the solar system planet's
   orbits. It takes the first of Kepler's laws of planetary
   motion which calculates the orbit's radius in a given moment
   for a given angle. If you don't know much of the subject
   I suggest you read and calculate wikipedia's page on the
   subject: http://en.wikipedia.org/wiki/Kepler's_laws_of_planetary_motion """


# TODO:
#
# - An abstract way of putting bodies in scene
# - Clean that mess down there (aka find a pythonic and clean way to do the job)
# - Planet sattelite support
# - A decent stat shower (in rt)
# - Clean the code up
# - User interaction
# - Maybe a shell and stuff
# - And lots of other stuff


######## MODULES ########

import visual as v
from math import *
from time import sleep
from thread import start_new_thread as nthread

#########################

####### CONSTANTS #######

G = 6.67e-11 # Univertial Gravitational Constant
pi = 3.141592653589793238462643 # pi | I could import it from math, but ...

#########################

####### SCENE CFG #######

v.scene.title = "Solar System"
v.scene.width = 1200
v.scene.height = 600

#########################

#def shell():
#   global planets, camera
#   
#   while True:
#      cmd = raw_input()
#      try:
#         cmd, target = cmd.split()
#      except:
#         pass
#   
#      if cmd == "center":
#         for planet in planets:
#            if planet == target:
#               camera = planet
#         if target == "sun":
#            camera = target
      

if __name__ == "__main__":
   
   planets = ["mercury", "venus", "earth", "mars", "jupiter"]
   
   
   """ data is a dictionary of lists of the following form [(semi-major axis, 0, 0), radius, color, aphelion, eccentricity, period (in days)]
      0 = pos
      1 = radius
      2 = color
      3 = aphelion
      4 = eccentricity
      5 = period
   """
   
   data = {"mercury":
      {"sm-axis": 57909100, "radius": 2439.7, "color": v.color.red, "aphelion": 69816900, \
      "eccentricity": 0.205630, "tperiod": 87.969, "rperiod": 58.646},
      
         "venus":
      {"sm-axis": 108208930, "radius": 6501.8, "color": v.color.white, "aphelion": 108942109, \
      "eccentricity": 0.0068, "tperiod": 224.70069, "rperiod": 243.0185},
      
         "earth":
      {"sm-axis": 147098074, "radius": 6371, "color": v.color.blue, "aphelion": 149597887.5, \
      "eccentricity": 0.016710219, "tperiod": 365.256366, "rperiod": 0.99726968 },
      
         "mars":
      {"sm-axis": 227939100, "radius": 3396, "color": v.color.orange, "aphelion": 249209300, \
      "eccentricity": 0.093315, "tperiod": 686.971, "rperiod": 1.025957},
      
         "jupiter":
      {"sm-axis": 778547200, "radius": 71492 / 10., "color": v.color.yellow, "aphelion": 816520800, \
      "eccentricity": 0.048775, "tperiod": 4311.572, "rperiod": 9.925 / 24} }

#         "moon":
#      { "sm-axis": 384399, "radius": 1737, "color": (0.5,0.5,0.5) , "aphelion": 405696, "eccentricity": 0.054, "tperiod": 27.321582 } }


   """ ps is a list with the values of p which is the semi-lactus rectum (what a weird name) of the ellipses """
   
   ps = {}
   for planet in data.keys():
      ps[planet] = data[planet]["aphelion"] * (1 - data[planet]["eccentricity"]**2)

      
   """ bodies is a list with vpython's spheres (the planets) """
   
   bodies = {}   
   for planet in data.keys():
      bodies[planet] = v.sphere(pos = (data[planet]["aphelion"],0,0) , radius = data[planet]["radius"] * 500, color = data[planet]["color"])


   """ trails is a list with vpython's curves (trails of the planets) """

   trails = {}
   for planet in data.keys():
#      trails[planet] = v.curve(color = data[planet]["color"])
      trails[planet] = v.curve(color = (1,1,1))
      
   """ tetas is a list with the values of teta which is the angle between the aphelion and r, the radius vector """
      
   tetas = {}
   alphas = {}
   for planet in data.keys():
      alphas[planet] = 0
      tetas[planet] = 0
   
   
   """ sun sphere object from vpython """
   sun = v.sphere(pos = (0,0,0), radius = 109 * 6371 * 10, color = v.color.red)
   
   
   """ xx and yy axis """
   xx = v.curve( x = (-800000000, 800000000) )
   yy = v.curve( y = (-800000000, 800000000) )
   v.label(pos = (800000000,0,0), text = "x")
   v.label(pos = (0,800000000,0), text = "y")
      
##   """ a new frame with the data of each planet """
##   data_frame = v.display(title = "Stats", width = 600, height = 200)
##   data_text = ""   
##
##   ############## LAME NOT SO PYTHONIC CODE ##############
##   # When I find a pythonic way to do it, I'll change it #
##   #######################################################
##
##   words = ["planet", "sm-axis", "radius", "aphelion", "eccentricity", "period (in days)"]
##
##   datas = [planets, [data[k][0][0] for k in data.keys()]]
##   for i in [1,3,4,5]:
##      datas.append( [ data[k][i] for k in data.keys() ] )
##   
##   for i in range(len(words)):
##      data_text += words[i] + "     "
##   data_text += "\n"
##      
##   for i in range(len(datas[0])):
##      for j in range(len(datas)):
##         data_text += str(datas[j][i]) + "      "
##      data_text += "\n"
##   
##   #######################################################
##   #######################################################
##   
##   v.label( display = data_frame, text = data_text )

   # Funny :p
   bodies["earth"].color = (1,1,1)
   bodies["earth"].material = v.materials.earth
   bodies["earth"].rotate(angle = pi/2, axis = (1,0,0))
   
   while True:
      v.rate(25)
      v.scene.center = bodies["earth"].pos
      
      for k in data.keys():
      
         bodies[k].rotate(angle = 2*pi / data[k]["rperiod"], axis = (0,0,1))#, axis = (0,1,0), origin = bodies[k].pos)
         tetas[k] += (2*pi) / data[k]["tperiod"] # add angular displacement
         
         r = ps[k] / (1 + data[k]["eccentricity"] * cos(tetas[k])) # radius
         

         bodies[k].x = cos(tetas[k]) * r
         bodies[k].y = sin(tetas[k]) * r
         
#         if tetas[k] >= 2*pi: # clean last orbit
#            tetas[k] = tetas[k] % (2*pi)
#            trails[k].pos = [(data[k]["aphelion"],0,0)]
         
         trails[k].append(pos = bodies[k].pos) # orbit's trail