MBA: Gás simples: mudanças entre as edições

De Física Computacional
Ir para navegação Ir para pesquisar
Sem resumo de edição
Sem resumo de edição
 
(Uma revisão intermediária pelo mesmo usuário não está sendo mostrada)
Linha 57: Linha 57:
= Distribuição de velocidades para um gás a certa temperatura =
= Distribuição de velocidades para um gás a certa temperatura =


A distribuição das velocidads das molécula de um é dada pela distribuição de Maxwell-Boltmann <ref>[https://edisciplinas.usp.br/pluginfile.php/48089/course/section/16461/qsp_chapter4-kineticTheory.pdf The Kinetic Theory of Gases]</ref>
A distribuição das velocidades das moléculas de um gás é dada pela distribuição de Maxwell-Boltzmann <ref>[https://edisciplinas.usp.br/pluginfile.php/48089/course/section/16461/qsp_chapter4-kineticTheory.pdf The Kinetic Theory of Gases]</ref>


<math display="block"> f\left(v\right)=\left(\frac{m}{2\pi kT}\right)^{3/2}4\pi v^{2}\exp\left(-\frac{mv^{2}}{kT}\right) </math>
<math display="block"> f\left(v\right)=\left(\frac{m}{2\pi kT}\right)^{3/2}4\pi v^{2}\exp\left(-\frac{mv^{2}}{kT}\right) </math>

Edição atual tal como às 20h54min de 9 de agosto de 2022

Anterior: Por que usar e o que são modelos baseados em indivíduos | Índice: Ecologia | Próximo: MBA: Caminhante aleatório

Colisão entre partículas

Considerando um modelo simples de gás ideal, não há força atuando sob as partículas, então a interação que ocorre entre as partículas se dá apenas por meio de colisões. Assim é necessário calcular a variação na velocidade de cada partícula após a colisão. Começando em uma dimensão, precisamos garantir a conservação do momento:

E da energia cinética:
Colocando o referencial de forma que , então as velocidades no novo referencial podem ser escritas como , de forma que da conservação do momento ficamos com:
Elevando ao quadrado:
Substituindo na equação de conservação de energia, podemos encontrar :
Calculando as raízes do segundo grau, temos:
E como uma solução é , mas queremos a situação em que , logo:
Substituindo em (1), temos:
Ou seja:
Retornando ao referencial original, sendo , temos então para :
E de maneira análoga para :
Um caso especial ocorre se então temos simplesmente:
Em duas dimensões, podemos reduzir o problema a uma dimensão, considerando que toda a alteração na velocidade devido a colisão entre partículas ocorre apenas na componente paralela a reta que liga o centro das duas esferas. Considerando que a posição de cada partícula é dada por , então um vetor entre as partículas pode ser escrito como . Podemos projetar ambas as velocidades então fazendo:
Obtemos o módulo da velocidade da partícula na direção e podemos trabalhar em uma única dimensão para encontrarmos a velocidade de ambas partículas após a colisão nesta dimensão. Ao fim podemos decompor novamente a velocidade final em ambos os eixos mantendo a mesma direção utilizando , onde é a diferença entre a posição das partículas e na componente .

Além disso vale lembrar que há a componente da velocidade perpendicular a , que vamos denotar como . Esta componente perpendicular permanece inalterada e pode ser visualizada na figura ao lado.

Colisão entre duas partículas mostrando explicitamente os vetores relacionados à partícula '"`UNIQ--postMath-00000023-QINU`"'.
Colisão entre duas partículas mostrando explicitamente os vetores relacionados à partícula 1.

Sendo assim, a velocidade final é dada por:

Onde é um vetor unitário que nos dá a direção entre os centros das partículas. Utilizando ad identidades:
ficamos então com:
Logo:
E uma vez que , então . Logo:
Ou ainda mais explícito, se fizermos , sendo as partículas e , onde , usando (2):
Temos então que:
Todo o cálculo exibido foi para uma partícula, para a segunda partícula, o cálculo é análogo.

Distribuição de velocidades para um gás a certa temperatura

A distribuição das velocidades das moléculas de um gás é dada pela distribuição de Maxwell-Boltzmann [1]


Foi realizado uma simulação com 100 partículas por 500.000 de passos em uma grade 100x100 com fronteira toroidal, o resultado encontra-se na figura abaixo assim como o código.

Resultado da simulação com a distribuição de Maxwell-Boltmann para kT=0.7 e m=1.


Código

Aqui chamamos a atenção que a principal diferença do código discutido para autômatos celulares é que agora não temos mais uma grade, mas um modelo contínuo. Também temos a possibilidade de configuração se as fronteiras serão toroidais (extremidades conectadas) ou não. Quanto aos resultados, precisamos considerar que é um modelo simplificado. O Calculo foi realizado para uma colisão entre duas partículas se dentro do mesmo passo deveria ocorrer mais de uma colisão, ou com mais de uma partícula simultaneamente o código falha, para remediar isso, o ideal é utilizar passos de tempo pequenos e/ou uma caixa relativamente grande em relação a quantidade de partículas. Além disso ao optarmos por uma fronteira toroidal ou não, estamos considerando uma caixa fechada ou não, ao considerarmos precisamos considerar a colisão com as fronteira também. Nesse último caso, uma simplificação simples é considerar que esta colisão ocorre primeiro. Novamente o problema surge se houver duas colisões no mesmo passo (agora podendo colidir com a fronteira), logo a recomendação é a mesma.

#Bibliotecas necessárias
from mesa import Agent, Model                                      #Classes Agente e Modelo
from mesa.time import SimultaneousActivation                       #Agendador simultâneo
from mesa.space import ContinuousSpace                             #Malha multigrid
import random                                                      #Número aleatórios
from mesa.datacollection import DataCollector                      #Importamos um coletor de dados
import matplotlib.pyplot as plt                                    #Plotar gráicos             
#from IPython.display import clear_output                           #Configurações de saída do Colab         
import numpy as np    

#AGENTE---------------------------------------------------------------------------------
class Agente(Agent):
    """Classe do agente"""
    def __init__(self, modelo,largura,altura):
        """Bibliotecas necessárias"""
        #modelo     - Modelo que ao qual o agente pertence
        super().__init__(self, modelo)        #Necessário para funcionar o modelo
        self.ppos     = (0,0)                 #Poição que irá se mover
        self.rai= 0.2                         #Raio
        self.lim=(largura,altura)             #Limites da caixa
        v=1                                   #Módulo da velocidade inicial
        ang = 2*random.random()*np.pi
        self.vel      = (v*np.cos(ang),v*np.sin(ang)) #Velocidade decomposta
        self.colisao=False
        self.ac=False
                
    def prox_pos(self):
      dt=0.01
      nx=self.pos[0]+self.vel[0]*0.01
      ny=self.pos[1]+self.vel[1]*0.01
      #Testar os limites
      if (nx < 0):
        nx = - nx; self.vel=(-self.vel[0],self.vel[1])
      if (nx > self.lim[0]):
        nx = 2*self.lim[0] - nx; self.vel=(-self.vel[0],self.vel[1])
      if (ny < 0):
        ny = - ny; self.vel=(self.vel[0], - self.vel[1])
      if (ny > self.lim[1]):
        ny = 2*self.lim[1] - ny;self.vel=(self.vel[0], - self.vel[1])

      return((nx,ny))
    
    def prox_vel(self):
        """Movimentação dos animais"""
        #Checamos se há colisão com outro animal:
        agentes = self.model.continuous.get_neighbors(pos=self.ppos, radius = 3*self.rai, include_center=True)
        for a in agentes:                                         #Percorremos a lista
          if (a.unique_id!=self.unique_id): 
              d=np.sqrt((a.ppos[0]-self.ppos[0])**2+(a.ppos[1]-self.ppos[1])**2)
              if (d<=2*self.rai):
                #Direção da colisão
                d=(a.ppos[0]-self.ppos[0],a.ppos[1]-self.ppos[1])/np.sqrt((a.ppos[0]-self.ppos[0])**2+(a.ppos[1]-self.ppos[1])**2)
                #Partícula 1
                u1=(a.vel[0]-self.vel[0])*d[0]+(a.vel[1]-self.vel[1])*d[1]
                nvel1=(self.vel[0]+u1*d[0],self.vel[1]+u1*d[1])
                #partícula 2,pela conservação de momento:
                nvel2=(self.vel[0]+a.vel[0]-nvel1[0],self.vel[1]+a.vel[1]-nvel1[1])   
                #Se nenhuma partícula tinha se colidido
                if (self.colisao==False and a.colisao==False):          
                  #Atualiza as velocidades
                  self.vel=nvel1; a.vel=nvel2
                  self.ac=True;a.ac=True   
                #Registra a colisão
                self.colisao=True; a.colisao=True       
                break
        return ()                                          

    def movimento(self):
      if(self.colisao):
        return (self.pos)
      else:
        return(self.ppos)
    
    def step(self):
        """Método obrigatório que prepara as mudanças"""
        self.ppos = self.prox_pos()
             
    def advance(self):
        """Método obrigatório que aplica as mudanças"""
        self.prox_vel()
        self.model.continuous.move_agent(self, self.movimento())

#MODELO
class Modelo(Model):
    """Modelo geral"""
    def gerar_particulas(self,N,largura,altura):
      """Função para gerar animais iniciais"""    
      #N        - Número de animais
      #largura  - Largura da grade
      #altura   - Altura da grade
      for n in range(N):
          con=True
          while(con):
            con=False
            a = Agente(self,largura,altura)
            X= largura*random.random()
            Y= altura*random.random()
            if (len(self.schedule.agents)==0):
              self.schedule.add(a)
              self.continuous.place_agent(a, (X, Y))
            else:
              agentes = self.continuous.get_neighbors(pos=(X,Y), radius = 2*a.rai, include_center=True)
              for viz in agentes:
                con = (True) if (len(agentes)==0) else (con)  #Precisamos checar se nasceu no mesmo lugar de outro animal
              if(con==False):
                self.schedule.add(a)
                self.continuous.place_agent(a, (X, Y))
      return()

    def __init__(self, modelo,N,seed=None):
        """Função chamada quando o modelo é inicializazdo"""
        # Modelo   - Dicionário com especificações do modelo
        # N        - Número de partículas
        # seed     - Seed dos números aleatórios do modelo do mesa
        
        # Desempacotar o dicionário
        largura = modelo["Largura"];altura=modelo["Altura"]; random.seed(modelo["Seed"])                     
        self.continuous = ContinuousSpace(x_max= largura, y_max=altura, torus= True, x_min=0, y_min= 0) #Configura a grade
        self.schedule   = SimultaneousActivation(self)            #Configura o agendador
        self.running   = True                                     #Condiçao para seguir executando o modelo
        self.gerar_particulas(N,largura,altura)    #Distribuir as ovelhas
           
    def step(self):
        """Avançar um passo do modelo"""
        self.schedule.step()                              #Avançamos os agentes
        #Limpamos as colisões
        for a in self.schedule.agents:                    #Vamos percorrer os agentes
          a.colisao=False
#Parâmetros
MAX =1000000
N=100
modelo  =  {"Largura":100     ,"Altura":100      ,"Seed":1} 
M = Modelo(modelo,N)
for i in range(MAX):
    M.step()
    if ((i+1)%(MAX/100)==0):
      #clear_output() 
      print(str(100*(1+i)/MAX)+"%")

vel=[]
E=0
for a in M.schedule.agents:
    vel.append(np.sqrt(a.vel[0]**2+a.vel[1]**2))
    E+=a.vel[0]**2+a.vel[1]**2
print(E)
with open('velocidades0.2R100K1000.txt', 'w') as f:
    f.write(str(vel))

Para plotar o gráfico, podemos então usar simplesmente:

import numpy as np
import matplotlib.pyplot as plt

a,b,c=plt.hist(vel,8, density=True, facecolor='g', alpha=0.75)
m=1
K=0.7
x=np.arange(0,2.5, 0.00001)
y = (4/np.sqrt(np.pi))*((m/K)**(3/2))*(x**2)*np.exp(-m*(x**2)/K)
plt.plot(x,y)

Citações


Anterior: Por que usar e o que são modelos baseados em indivíduos | Índice: Ecologia | Próximo: MBA: Caminhante aleatório