Simulação-Secagem-PVC-Dados-INEOS-BAEYENS-Braskem-SA

import numpy as np from scipyintegrate import solveivp 1 CONSTANTES UNIVERSAIS E FÍSICAS Constantes físicas g 981 ms² Gravidade Rg 8314 JmolK Constante universal dos gases Mwair 0029 kgmol Massa molar do ar seco MwH2O 0018 kgmol Massa molar da água Propriedades gerais do PVC não específicas dos casos Cpp 1670 JkgK Calor específico do PVC seco lambdap 017 WmK Condutividade térmica do PVC CPS 1670 JkgK Calor específico da partícula redundante Conversão de temperatura Tconv 27315 K Para converter C em Kelvin Condições comuns TS0 30 Tconv K Temperatura inicial da partícula Up0 146 ms Velocidade inicial da partícula estimado Ug0 112 ms Velocidade inicial do gás estimado Pop0 101325 Pa Pressão atmosférica Twall 50 Tconv K Temperatura da parede assumida hwall 250 WmK Perda térmica por metro de secador 2 Propriedades Específicas do Caso INEOS Propriedades do PVC INEOS dpd 160e6 m Diâmetro médio das partículas rhop 1116 kgm³ Densidade do PVC chi 015 Porosidade tau 30 Tortuosidade XS0 033 kgkg Umidade inicial da partícula Geometria INEOS Dpipe 17 m Diâmetro interno do secador L 220 m Comprimento do secador Operação INEOS TG0 180 Tconv K Temperatura do ar de entrada YG0 0004 35 kgkg Umidade do ar de entrada Fp0 8545 1000 3600 kgs Vazão de PVC úmido Fg0 45 1000 3600 kgs Vazão de ar seco relação Fg0 Fp0 5266 3 Propriedades Específicas do Caso BAEYENS Propriedades do PVC Baeyens dpd 180e6 m Diâmetro médio das partículas rhop 1116 kgm³ Densidade do PVC chi 015 Porosidade tau 80 Tortuosidade XS0 026 kgkg Umidade inicial da partícula Geometria Baeyens Dpipe 125 m Diâmetro interno do secador L 250 m Comprimento do secador Operação Baeyens TG0 127 Tconv K Temperatura do ar de entrada YG0 0003 kgkg Umidade do ar de entrada Fp0 6670 1000 3600 kgs Vazão de PVC úmido Fg0 4648 1000 3600 kgs Vazão de ar seco Relação Fg0 Fp0 6968 3 Caso Braskem SA PVCAL02 Propriedades do PVC dpd 15e4 m Diâmetro médio das partículas rhop 812 kgm³ Densidade do PVC chi 035 Porosidade tau 30 Tortuosidade XS0 033 kgkg Umidade inicial da partícula Geometria Baeyens L 26934 m Comprimento do secador Operação Baeyens TG0 180 Tconv K Temperatura do ar de entrada YG0 0004 kgkg Umidade do ar de entrada Fp0 30 1000 3600 kgs Vazão de PVC seco Fg0 180 1000 3600 kgs Vazão de ar seco relação Fg0 Fp0 6 CPL 4181 JkgK calor específico da água líquida CPA 1012 JkgK INDEXADOR PARA CORRELAÇÃO DE NUSSELT UTILIZADA Observação As correlações disponíveis são 1 RanzMarshall LeviHevroni et al 1995 2 Baeyens et al 1995 3 De Brandt Baeyens et al 1995 4 Gamson Kerker and Terwiesch 1985 5 Khotari Kunii and Levenspiel 1969 PIOR SEMPRE 6 Bandrowski 1977 idx 6 2 FUNÇÕES AUXILIARES PROPRIEDADES DO GÁS def CalculateRaioz Retorna o raio interno em metros ao longo do secador com base na interpolação dos diâmetros definidos nos trechos geométricos do secador trechos 0000 2577 10 10 Diâmetro constante 10 m 2577 9310 10 20 Convergente 10 m 20 m 9310 16790 20 20 Diâmetro constante 20 m 16790 18809 20 17 Convergente 20 m 17 m 18809 26934 17 12 Convergente 17 m 12 m for zi zf di df in trechos if zi z zf if di df return di 2 diâmetro constante raio else L zf zi dinterp di df di z zi L return dinterp 2 converte para raio def CalculateRaioz return 125 2 def CalculateAcRatio Retorna a área da seção transversal do secador pneumático em função de z return nppi Ratio2 def CalculatedHvapTS Calcula o calor latente de vaporização da água dHvap com base em um ajuste polinomial de grau 3 Entrada TS float ou arraylike Temperatura da partícula em Kelvin K Saída dHvap float ou arraylike Calor latente em Jkg TC TS 27315 Conversão de Kelvin para Celsius Coeficientes do polinômio grau 3 ajuste baseado em dados tabelados a0 2501455505 a1 2179152 a2 0002803 a3 0000001 dHvapkJkg a0 a1TC a2TC2 a3TC3 return dHvapkJkg 1000 Convertendo para Jkg def CalculaterhoLTS Calcula a densidade da água líquida rhoL em kgm³ com base em um ajuste polinomial de grau 2 feito sobre dados experimentais no intervalo de 10 C a 110 C 28315 K a 38315 K Entrada TS float ou arraylike Temperatura em Kelvin K Saída rhoL float ou arraylike Densidade da água líquida em kgm³ TC TS 27315 Conversão para Celsius Coeficientes do polinômio de grau 2 ajuste com seus dados a 0003155 b 0121666 c 100186 return a TC2 b TC c kgm³ def CalculateCpgTG CONFERIDO Calor específico do ar seco JkgK ajustado por polinômio de grau 3 baseado em dados tabelados entre 50C e 210C Entrada TG float Temperatura do gás K Saída Cpg float Calor específico ajustado JkgK return 147306397e07 TG3 324804744e04 TG2 225536353e01 TG 959548912 def CalculateCpvTS CONFERIDO Calor específico do vapor de água JkgK ajustado por polinômio de grau 3 baseado em dados tabelados entre 50C e 210C Entrada TS float Temperatura da partícula úmida K Saída Cpv float Calor específico ajustado do vapor JkgK return 526094276e07 TS3 810256584e04 TS2 219725199e01 TS 186834356 def CalculatemugTG CONFERIDO Viscosidade dinâmica do ar seco Pas com ajuste polinomial grau 2 baseado em dados do NIST para o intervalo 323483 K Entrada TG float Temperatura do gás K Saída mug float Viscosidade dinâmica do ar Pas return 603354978e11 TG2 100731845e07 TG 729726940e06 def CalculatelambdagTG CONFERIDO Condutividade térmica do ar seco WmK com ajuste polinomial grau 2 baseado em dados reais no intervalo de 323483 K Entrada TG float Temperatura do gás K Saída lambdag float Condutividade térmica ajustada WmK return 197510823e08 TG2 741752976e05 TG 548126241e03 def CalculateDairH2OTG Difusividade do vapor de água no ar seco m²s usando a correlação de Fuller Entrada TG float ou arraylike Temperatura do gás em Kelvin K Saída Difusividade do vapor de água no ar m²s Fórmula D 11757 x 10 ⁹ x T175 Fonte Fuller et al return 11757e9 TG175 def CalculatePrCpg mug lambdag CONFERIDO Calcula o número de Prandtl do gás Entradas Cpg float Calor específico do gás JkgK mug float Viscosidade dinâmica Pas lambdag float Condutividade térmica WmK Saída Pr float Número de Prandtl Observação Relação entre difusão de quantidade de movimento e de calor return Cpg mug lambdag def CalculateScmug rhog Dab OK Calcula o número de Schmidt do gás Entradas mug float Viscosidade dinâmica Pas rhog float Densidade do gás kgm³ D float Difusividade do vapor no gás m²s Saída Sc float Número de Schmidt Observação Relação entre difusão de quantidade de movimento e de massa return mug rhog Dab def CalculaterhogTG YG Pop CONFERIDO Calcula a densidade do gás úmido como mistura de ar seco vapor de água Entradas TG float Temperatura do gás K YG float Umidade absoluta do gás kg vapor kg ar seco Pop float Pressão de operação Pa Saída rhog float Densidade da mistura gasosa úmida kgm³ Observação Baseado na equação dos gases ideais aplicada a misturas binárias Usa massas molares em kgmol para consistência com o restante do modelo MA MwH2O kgkmol vapor de água MB Mwair kgkmol ar seco denominador YG MA 1 MB numerador Pop Rg TG 1 YG return numerador denominador 3 FUNÇÃO PARA COEFICIENTE DE ARRASTE Wen Yu def CalculateRerhog Ur dp mug OK Calcula o número de Reynolds para uma partícula esférica em escoamento Entradas rhog float Densidade do gás kgm³ Ur float Velocidade relativa entre gás e partícula ms dp float Diâmetro da partícula m mug float Viscosidade dinâmica do gás Pas Saída Re float Número de Reynolds Observação Representa a razão entre forças inerciais e viscosas no escoamento ao redor da partícula return rhog Ur dp mug def CalculateCdRe alphap CONFERIDO Calcula o coeficiente de arraste Cd com base na correlação de Wen Yu incluindo o efeito da fração de fase gasosa alphag17 conforme a Equação 29 do artigo de Aubin et al 2013 Entradas Re float Número de Reynolds baseado na partícula alphap float Fração volumétrica de partículas Saída Cd float Coeficiente de arraste adimensional alphag 10 alphap fração volumétrica da fase gasosa Cdbase 24 Re 1 015 Re0687 return Cdbase alphag17 def CalculateRepiperhog Ug Dpipe mug CONFERIDO Calcula o número de Reynolds do escoamento do gás no tubo do secador Entradas rhog float Densidade do gás úmido kgm³ Ug float Velocidade do gás ms Dpipe float Diâmetro interno do tubo m mug float Viscosidade dinâmica do gás Pas Saída Repipe float Número de Reynolds no tubo return rhog Ug Dpipe mug def CalculatefwallRepipe Calcula o fator de atrito com a parede fwall com base no número de Reynolds do escoamento no tubo Repipe conforme as equações do artigo de Aubin et al 2013 usando interpolação linear na faixa de transição Entrada Repipe float Número de Reynolds do gás no tubo do secador Saída fwall float Fator de atrito adimensional Regime laminar if Repipe 2100 return 64 Repipe Regime turbulento elif Repipe 10000 return 00791 Repipe025 Faixa de transição interpolação linear else flaminar 64 2100 fturbulento 00791 10000025 Interpolação fwall flaminar fturbulento flaminar Repipe 2100 10000 2100 return fwall def CalculatefwallRepipe Calcula o fator de atrito usando a correlação de Churchill 1977 válida para todos os regimes laminar transição turbulento A 2457 nplog1 7Repipe09 02700001Dpipe16 B 37530Repipe16 fwall 8 8Repipe12 1A B15112 return fwall 4 FUNÇÕES AUXILIARES TRANSFERÊNCIA E SECAGEM def calculatedpXS dpd rhop rhoL Xcrit CONFERIDO Calcula o diâmetro da partícula úmida com base na condição fornecida Args XS float Umidade da partícula escalar dpd float Diâmetro da partícula seca rhop float Densidade da partícula seca rhoL float Densidade da água líquida Xcrit float Umidade crítica Returns float Diâmetro úmido da partícula if XS Xcrit return dpd 10 rhop rhoL Xcrit rhop rhoLXS 1 3 else return dpd def CalculatealphapFp0 rhop Ac Up ADIÇÃO EXTERNA AO ARTIGO Calcula a fração volumétrica de partículas alphap dinamicamente com base na vazão mássica velocidade e geometria do tubo Entradas Fp0 float Vazão mássica de partículas secas kgs rhop float Densidade das partículas kgm³ Ac float Área da seção transversal do tubo m² Up float Velocidade atual das partículas ms Saída alphap float Fração volumétrica de partículas no tubo return Fp0 rhop Ac Up def CalculateNuidx Re Pr alphap Calcula o número de Nusselt Nu usando diferentes correlações conforme o índice idx fornecido Entradas idx int Índice que determina qual correlação será usada 1 a 6 Re float Número de Reynolds da partícula Pr float Número de Prandtl alphap float Porosidade da partícula Saída Nu float Número de Nusselt adimensional Observação As correlações disponíveis são 1 RanzMarshall LeviHevroni et al 1995 2 Baeyens et al 1995 3 De Brandt Baeyens et al 1995 4 Gamson Kerker and Terwiesch 1985 5 Khotari Kunii and Levenspiel 1969 6 Bandrowski 1977 if idx 1 Nu 2 06 Re05 Pr13 elif idx 2 Nu 015 Re elif idx 3 Nu 016 Re13 Pr067 elif idx 4 Nu 106 Re059 Pr13 elif idx 5 Nu 0003 Re13 elif idx 6 Nu 000114 alphap05984 Re08159 else raise ValueErrorÍndice inválido Use um valor de 1 a 6 para selecionar a correlação return Nu def CalculateCphCpp XS Cpv CONFERIDO Calcula o calor específico da partícula úmida Cph conforme o artigo Entradas Cpp float Calor específico do PVC seco JkgK XS float Umidade da partícula kg água kg sólido seco Cpv float Calor específico da água fase vapor JkgK Saída Cph float Calor específico da partícula úmida JkgK return Cpp XS Cpv def CalculatehNu lambdag dp OK Calcula o coeficiente de transferência de calor convectivo h Entradas Nu float Número de Nusselt lambdag float Condutividade térmica do gás WmK dp float Diâmetro da partícula m Saída h float Coeficiente de transferência de calor Wm²K return Nu lambdag dp def Calculatekyxh Cpg Sc Pr CONFERIDO Calcula o coeficiente de transferência de massa ky usando a analogia de ChiltonColburn Entradas h float Coeficiente de transferência de calor Wm²K Cpg float Calor específico do gás JkgK Sc float Número de Schmidt Pr float Número de Prandtl Saída kyx float Coeficiente de transferência de massa kgm²s return h Cpg Sc Pr23 def CalculateDappDab chi tau CONFERIDO Calcula a difusividade aparente no interior da partícula Dapp usando a correlação do artigo de Aubin et al 2013 Entradas Dab float Difusividade do vapor de água no ar m²s chi float Porosidade da partícula tau float Tortuosidade interna da partícula Saída Dapp float Difusividade aparente m²s return Dab chi tau def CalculateBiMkyx dp Dapp rhog OK Calcula o número de Biot para transferência de massa BiM conforme a equação 7 do artigo de Aubin et al 2013 Entradas ky float Coeficiente de transferência de massa kgm²s dp float Diâmetro da partícula m Dapp float Difusividade aparente no interior da partícula m²s rhog float Densidade do gás úmido kgm³ Saída BiM float Número de Biot para transferência de massa R dp 2 return kyx R Dapp rhog def CalculateXcritrhowater chi rhop CONFERIDO Calcula a umidade crítica Xcrit da partícula Entradas rhowater float Densidade da água líquida kgm³ chi float Porosidade da partícula rhop float Densidade da partícula kgm³ Saída Xcrit float Umidade crítica kg água kg sólido seco Observação Representa o limite abaixo do qual a difusão interna controla a secagem rhowater chi rhop return rhowater chi rhop def CalculateRhowetrhops XS dpd dp Calcula a densidade específica da partícula úmida rhodry corrigida pelo inchaço ou retração da partícula durante a secagem ou hidratação Entradas rhops float Densidade da partícula seca kgm³ XS float Umidade da partícula kg água kg sólido seco dpd float Diâmetro da partícula seca m dp float Diâmetro da partícula úmida m Saída rhodry float Densidade específica da partícula úmida kgm³ Observação A equação assume que a massa total da partícula é aumentada pela fração de água presente 1 XS e que o volume da partícula muda com o cubo da razão de diâmetros dpd dp³ return rhops 1 XS dpd dp 3 def CalculateKyXS XS0 kyx BiM rhop rhowater chi Xcrit Calcula o coeficiente global de transferência de massa Ky Entradas XS kgkg Umidade da partícula XS0 kgkg Umidade inicial da partícula kyx ms Coeficiente de transferência externa de massa BiM Número de Biot modificado para massa rhop kgm³ Densidade da partícula rhowater kgm³ Densidade da água líquida chi Porosidade da partícula Xcrit kgkg Umidade crítica limite para mudança de regime difusional Saída Ky ms Coeficiente global de transferência de massa Regime difusional interno significativo if XS Xcrit return kyx BiM XS XS013 1 1 1 XS0 Xcrit Regime com controle difusional reduzido else termo rhop XS rhowater 1 chi return kyx termo 2 3 def CalculatePwsTS Calcula a pressão de vapor saturante da água pws usando o ajuste Antoine modificado Entrada TS temperatura da partícula K Saída pws pressão de vapor saturante Pa TC TS 27315 conversão para C return npexp1683 TC 22874 TC 64075 def CalculateYsatPopPws Calcula a umidade absoluta de saturação Ysat da partícula a uma dada temperatura TS K Entradas TS temperatura da partícula K Pop pressão de operação do sistema Pa Saída Ysat umidade absoluta de saturação kg vapor kg ar seco return 0622 Pws Pop Pws def CalculateIgkalphap rhog Ur Cd dp Calcula a interação de quantidade de movimento do gás para a partícula Igp Entradas alphap float Fração volumétrica de partículas rhog float Densidade do gás kgm³ Ur float Velocidade relativa gás partícula ms Cd float Coeficiente de arraste dp float Diâmetro da partícula m Saída Igp float Interação de quantidade de movimento do gás para a partícula kgm²s² return 3 alphap rhog Ur2 Cd 4 dp def Calculatefwgfwall rhog Ug Dpipe Calcula a força de atrito da parede sobre o gás Ffwg Entradas f float Fator de atrito adimensional rhog float Densidade do gás kgm³ Ug float Velocidade do gás ms Dpipe float Diâmetro interno do tubo m Saída fwg float Força de atrito da parede sobre o gás kgm²s² return fwall rhog Ug2 2 Dpipe def CalculateNwKy Ysat YG Calcula o fluxo de massa evaporada Nw da partícula para o gás Entradas Ky float Coeficiente global de transferência de massa kgm²s Ysat float Umidade de saturação do gás kgkg YG float Umidade absoluta atual do gás kgkg Saída Nw float Fluxo de massa evaporada kgm²s return Ky Ysat YG def CalculateGammapalphap sp Nw Calcula a fonte de massa volumétrica Gamma no gás devido à evaporação da partícula Entradas alphap float Fração volumétrica de partículas sp float Área superficial específica da partícula 1m Nw float Fluxo de massa evaporada kgm²s Saída Gamma float Fonte de massa volumétrica kgm³s return alphap sp Nw def CalculateGammagGammap Calcula a fonte de massa volumétrica Gamma no gás devido à evaporação da partícula Entradas alphap float Fração volumétrica de partículas sp float Área superficial específica da partícula 1m Nw float Fluxo de massa evaporada kgm²s Saída Gamma float Fonte de massa volumétrica kgm³s return Gammap def calculateHgwdHvap TG CpvTS Calcula a entalpia da água evaporada calor de vaporização mais aquecimento do vapor até a temperatura do gás onde ele se mistura Entradas dHvap float Calor latente de vaporização da água Jkg TG float Temperatura do gás K TS float Temperatura da partícula K Cpv float Calor específico do vapor de água JkgK Saída Hgw float Entalpia da água evaporada Jkg return dHvap Cpv TG TS def calculateQgtopalphap Ac sp h TG TS Calcula o fluxo de calor do gás para a partícula Entradas alphap float Fração volumétrica de partículas Ac float Área da seção transversal m² sp float Área superficial específica da partícula 1m h float Coeficiente de transferência de calor Wm²K TG float Temperatura do gás K TS float Temperatura da partícula K Saída Qgtop float Fluxo de calor do gás para a partícula Wm return alphap Ac sp h TG TS def calculateQptogQgtop Calcula o fluxo de calor da partícula para o gás Entrada Qgtop float Fluxo de calor do gás para a partícula Wm Saída Qptog float Fluxo de calor da partícula para o gás Wm sinal oposto de Qgtop return Qgtop def calculateQwtoghwall Twall TG Calcula o fluxo de calor da parede para o gás Entradas hwall float Coeficiente de troca térmica com a parede Wm²K Twall float Temperatura da parede K TG float Temperatura do gás K Saída Qwtog float Fluxo de calor da parede para o gás Wm return hwall Twall TG zprintproximo 00 Começa no 0 m passoprint 10 Intervalo entre prints 1 metro 5 FUNÇÃO COM O SISTEMA DE EDOs MODELO COMPLETO def dryermodelz y params global zprintproximo Variáveis de estado XS YG TS TG Up Ug Pop y XS0 Fg0Fp0Twallrhop params Ratio CalculateRaioz Ac CalculateAcRatio Dpipe Ratio 2 1 Propriedades termofísicas calculadas a partir da temperatura do gáspartícula Cpg CalculateCpgTG Cpv CalculateCpvTS mug CalculatemugTG lambdag CalculatelambdagTG dHvap CalculatedHvapTS Dab CalculateDairH2OTG rhog CalculaterhogTG YG Pop rhoL CalculaterhoLTS Xcrit CalculateXcritrhoL chi rhop dp calculatedpXS dpd rhop rhoL Xcrit alphap CalculatealphapFp0 rhop Ac Up alphag 1 alphap sp 6 dp 2 Propriedades adimensionais e coeficientes de transferência Pr CalculatePrCpg mug lambdag Sc CalculateScmug rhog Dab Ur Ug Up Re CalculateRerhog absUr dp mug Nu CalculateNuidx Re Pr alphap h CalculatehNu lambdag dp kyx Calculatekyxh Cpg Sc Pr Dapp CalculateDappDab chi tau BiM CalculateBiMkyx dp Dapp rhog Ky CalculateKyXS XS0 kyx BiM rhop rhoL chi Xcrit Pws CalculatePwsTS Ysat CalculateYsatPopPws Cd CalculateCdRe alphap Igp CalculateIgkalphap rhog Ur Cd dp Ipg Igp Repipe CalculateRepiperhog Ug Dpipe mug fwall CalculatefwallRepipe fwg Calculatefwgfwall rhog Ug Dpipe 3 Fluxo de evaporação e fonte de massa Nw CalculateNwKy Ysat YG Gammap CalculateGammapalphap sp Nw Gammag CalculateGammagGammap Hgw calculateHgwdHvap TG CpvTS Qgtop calculateQgtopalphap Ac sp h TG TS Qptog Qgtop Qwtog calculateQwtoghwall Twall TG 5 CÁLCULOS DOS BALANÇOS MASSA ENERGIA E MOMENTO BALANÇOS DE MASSA dXSdz Ac Gammap Fp0 negativo umidade decresce dYGdz Ac Gammag Fg0 BALANÇO DE ENERGIA dTSdz Qgtop Ac Gammap Hgw Fp0 CPS XS CPL dTGdz Qptog Qwtog Ac Gammag Hgw Fg0 Cpg YG Cpv BALANÇO DE ENERGIA BALANÇO DE MOMENTO Derivadas internas aproximadas ou calculadas previamente dalphapdz 0 ou cálculo baseado no Up se desejar depois dalphagdz dalphapdz drhopdz 0 PVC drhogdz 0 pode implementar o cálculo exato depois Construção da matriz A e vetor B A nparray 2 alphag rhog Ug 0 alphag 0 2 alphap rhop Up alphap alphag rhog alphap rhop 0 B nparray dalphagdz rhog Ug2 alphag drhogdz Ug2 alphag rhog g Ipg Up Ug Gammag fwg dalphapdz rhop Up2 alphap drhopdz Up2 alphap rhop g Igp rhog Ug dalphagdz alphag Ug drhogdz rhop Up dalphapdz alphap Up drhopdz Resolver o sistema derivatives nplinalgsolveA B dUgdz dUpdz dPdz derivatives DIAGNOSTIOOS if z zprintproximo zprintproximo 05 printfAc Gammap Hgw Ac Gammap Hgw printfQgtop Qgtop if XS Xcrit printfXS Xcrit DIAGNOSTIOOS return dXSdz dYGdz dTSdz dTGdz dUpdz dUgdz dPdz 6 CONDIÇÕES INICIAIS E INTEGRAÇÃO DO MODELO Parâmetros do modelo passados para dryermodel params XS0 Fg0Fp0Twallrhop Vetor de variáveis de estado iniciais y0 XS0 YG0 TS0 TG0 Up0 Ug0 Pop0 sol solveivp fundryermodel Sistema de EDOs tspan0 L Intervalo de integração z 0 até z L y0y0 Condições iniciais argsparams Parâmetros adicionais denseoutputTrue Permite interpolação contínua com solsolz methodLSODA Integrador escolhido rtol1e5 atol1e6 maxstep001 passo máximo de 10 cm ultimovalorxs soly01 ultimovalorts soly21 ultimovalortg soly31 printfÚmidade final ultimovalorxs printfTemperatura Final G ultimovalortg 27315 printfTemperatura Final S ultimovalorts 27315 import matplotlibpyplot as plt Eixo espacial para plotagem suave zplot nplinspace0 L 200 XS YG TS TG Up Ug Pop solsolzplot Inclui todas as variáveis Converte temperaturas para Celsius TSC TS 27315 TGC TG 27315 Converte úmidades Converte a pressão de Pa para atm Popatm Pop 101325 Calcula a velocidade relativa Ug Up Ur Ug Up Cria figura com 2 linhas e 3 colunas fig axs pltsubplots2 3 figsize18 8 Gráfico 1 Umidade da partícula XS axs0 0plotzplot XS colorblue axs0 0settitleUmidade da partícula XS axs0 0setxlabelComprimento z m axs0 0setylabelXS kgkg axs0 0gridTrue Gráfico 4 Temperatura da partícula TS em C axs1 0plotzplot TSC colorred axs1 0settitleTemperatura da partícula TS axs1 0setxlabelComprimento z m axs1 0setylabelTS C axs1 0gridTrue Gráfico 2 Umidade do gás YG axs0 1plotzplot YG colorgreen axs0 1settitleUmidade do gás YG axs0 1setxlabelComprimento z m axs0 1setylabelYG kgkg axs0 1gridTrue Gráfico 3 Temperatura do gás TG em C axs0 2plotzplot TGC colororange axs0 2settitleTemperatura do gás TG axs0 2setxlabelComprimento z m axs0 2setylabelTG C axs0 2gridTrue Gráfico 5 Pressão do gás P em atm axs1 1plotzplot Popatm colorpurple axs1 1settitlePressão do gás P axs1 1setxlabelComprimento z m axs1 1setylabelPressão atm axs1 1gridTrue Gráfico 6 Velocidade Relativa Ug Up axs1 2plotzplot Ur colorblack axs1 2settitleVelocidade Relativa Ug Up axs1 2setxlabelComprimento z m axs1 2setylabelUr ms axs1 2gridTrue Ajustes finais for ax in axsflat axticklabelformatstyleplain useOffsetFalse plttightlayout pltshow