Optimización del código del simulador de reactores de fusión nuclear Pretor

Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 1 Resumen Este proyecto consiste en el desarrollo de una versión preliminar de un simulador de fusión nuclear en entorno Windows llamado WinPretor. Esta aplicación informática se basa en el código fuente del programa PRETOR PC desarrollado en MSDOS y lenguaje Basic. La aplicación PRETOR-PC adolece de muchas carencias, sobre todo en el aspecto gráfico, dada su naturaleza MS DOS. Con la versión WinPretor, se pretende mejorar y ampliar las capacidades gráficas aprovechando al máximo las posibilidades que nos brinda el entorno Windows. Para desarrollar las bases de la aplicación WinPretor, primero ha sido necesario entender el funcionamiento de la aplicación MSDOS. Una tarea realmente complicada, ya que el código original está absolutamente indocumentado. Por lo tanto, se han documentado las partes más importantes del código antes de proceder a la nueva versión Windows. El proyecto tiene una primera parte introductoria donde se explica el fenómeno físico de la fusión nuclear; a continuación se comenta brevemente el trasfondo físico y matemático en el que se sustenta el código de Pretor PC. En los siguientes capítulos se analiza a fondo la estructura, diseño y funcionalidades del código fuente y se explican los pasos seguidos en el desarrollo de la nueva versión Windows con sus mejoras y nuevas funcionalidades. En el último capítulo se mencionan futuros proyectos que se podrían emprender utilizando éste como base. Pág. 2 Memoria Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 3 Sumario RESUMEN ___________________________________________________1 SUMARIO ____________________________________________________3 1. INTRODUCCIÓN __________________________________________7 1.1. Objetivos del proyecto......................................................................................7 1.2. Alcance del proyecto........................................................................................7 2. CONTEXTO ENERGÉTICO ACTUAL __________________________8 3. LA FUSIÓN NUCLEAR ____________________________________10 3.1. El fenómeno físico..........................................................................................10 3.2. Confinamiento magnético..............................................................................14 3.3. Reactores Tokamak.......................................................................................17 3.4. Balance de energía........................................................................................20 4. MODELO PRETOR________________________________________23 4.1. Estructura del Código PRETOR PC..............................................................25 4.1.1. Ficheros...............................................................................................................25 4.1.2. Funciones ............................................................................................................26 4.1.3. Acceso a ficheros................................................................................................27 4.1.4. Estructuras de datos............................................................................................30 4.1.5. Variables Globales significativas.........................................................................34 4.1.6. Diagrama de flujo simplificado del motor de simulación ....................................39 5. WINPRETOR ALFA _______________________________________43 5.1. Fases del desarrollo.......................................................................................44 5.2. Módulos funcionales ......................................................................................46 5.2.1. Módulo de entorno gráfico...................................................................................47 5.2.2. Módulo de carga de ficheros...............................................................................49 5.2.3. Módulo de simulación simplificado......................................................................53 5.2.4. Pantalla de visualización de parámetros............................................................54 5.2.5. Módulo de visualización de variables del plasma...............................................55 5.2.6. Módulo de perfiles de plasma.............................................................................56 5.2.7. Módulo de evoluciones temporales ....................................................................58 6. FUTUROS DESARROLLOS ________________________________60 6.1. Winpretor versión definitiva............................................................................60 6.1.1. Acceso completo a ficheros................................................................................60 Pág. 4 Memoria 6.1.2. Mejorar el control de la simulación.....................................................................60 6.1.3. Soporte completo para eventos..........................................................................60 6.1.4. Compatibilidad absoluta con Excel.....................................................................60 6.1.5. Optimizar el código de simulación......................................................................61 6.1.6. Implementar rutinas de cálculo externas al código tipo DLL .............................61 6.1.7. Desarrollo de una versión cliente-servidor.........................................................61 6.1.8. Desarrollo de una versión para computación distribuida...................................63 6.1.9. Desarrollo de una versión docente simplificada.................................................63 CONCLUSIONES _____________________________________________64 AGRADECIMIENTOS__________________________________________65 BIBLIOGRAFÍA_______________________________________________67 Referencias bibliográficas .......................................................................................67 Bibliografía complementaria:...................................................................................68 Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 5 ANEXOS ÍNDICE _______________________________________________________ A. PRESUPUESTO ECONÓMICO ________________________________ A.2. Recursos humanos............................................................................................ A.2.1. Análisis funcional de la aplicación Pretor PC: 25 h................................................ A.2.2. Análisis de la estructura del código fuente Pretor PC: 50h.................................... A.2.3. Documentación del código Pretor PC: 150h .......................................................... A.2.4. Formación sobre los fundamentos físicos subyacentes: 200h.............................. A.2.5. Separación de las rutinas de simulación y de interacción con el entorno gráfico y sistema operativo: 150 h ...................................................................................... A.2.6. Reconstrucción de rutinas y corrección de errores de sintaxis: 100h.................... A.2.7. Diseño conceptual de la aplicación WinPretor: 50 h.............................................. A.2.8. Adaptación del módulo de simulación Pretor PC a WinPretor: 100h .................... A.2.9. Desarrollo del entorno gráfico WinPretor: 120h...................................................... A.2.10. Testeo y depuración básica: 60h............................................................................ A.2.11. Tareas administrativas: 90 h................................................................................... A.3. Recursos materiales .......................................................................................... A.3.1. Ordenador Personal................................................................................................ A.3.2. Otros ........................................................................................................................ A.4. Cómputo total..................................................................................................... B. VARIABLES _______________________________________________ B.1. Matriz de variables PHYS (n,m)......................................................................... B.2. Matriz INF(n,m).................................................................................................. B.3. Vector Pval(n) y Ptitle(n).................................................................................... C. CÓDIGO FUENTE __________________________________________ C.1. Definción de variables globales......................................................................... C.2. Nuevas Clases................................................................................................... C.2.1. Clase OclsRead....................................................................................................... C.2.2. Clase OClsProject................................................................................................... C.2.3. Clase OclsGraph..................................................................................................... C.3. Nuevas Funciones ............................................................................................. C.3.1. Print2Excel............................................................................................................... C.3.2. LoadProject.............................................................................................................. C.3.3. Step3........................................................................................................................ Pág. 6 Memoria C.3.4. Step2 ....................................................................................................................... C.3.5. LoadInterfaceConfig................................................................................................ C.3.6. SaveInterfaceConfig................................................................................................ C.3.7. LoadDat con algoritmo de interpolación nodal....................................................... C.3.8. OmsgBox................................................................................................................. C.4. Rutinas de acceso a Datas modificadas ........................................................... C.4.1. InitDatas................................................................................................................... C.4.2. CargarPhysname.................................................................................................... C.4.3. CargarPanels y CargarInfSub................................................................................. C.4.4. CargarInf.................................................................................................................. D. FORMULARIOS MÁS IMPORTANTES __________________________ D.1. Formulario de Evoluciones temporales OfrmEVTemp...................................... D.2. Formulario de perfiles de plasma OfrmPerfil..................................................... D.3. Formulario OfrmPlot........................................................................................... D.4. Formulario OfrmBalance.................................................................................... D.5. Formulario ofrmOpenProject ............................................................................. Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 7 1. Introducción A continuación se exponen los objetivos y el alcance de este proyecto. 1.1. Objetivos del proyecto Los simuladores informáticos son una herramienta común en la fusión nuclear y en general, en todos los campos donde la realización empírica implica grandes costes. En la actualidad, existen multitud de simuladores de fusión nuclear, la mayoría de ellos adolecen de falta de flexibilidad y tienen entornos gráficos difíciles de utilizar y complicados sistemas para la actualización de sus rutinas de cálculo. El objetivo principal de este proyecto es desarrollar las bases de la aplicación WINPRETOR, un simulador que permita la máxima flexibilidad, un entorno gráfico intuitivo y que facilite las actualizaciones del código de cálculo. El sistema operativo utilizado será Windows XP y el código fuente se desarrollará en Visual Basic. Se partirá del código del programa PRETOR PC para desarrollar el nuevo simulador WINPRETOR. Un primer paso para llegar a este fin es entender la física subyacente y analizar en profundidad el código fuente PRETOR PC. El segundo paso es desarrollar la versión Windows modularmente desde la versión MSDOS. 1.2. Alcance del proyecto El proyecto pretende desarrollar los aspectos más básicos de WINPRETOR para que sirva como base para un futuro desarrollo completo. Por lo tanto, dado el carácter preliminar de la versión la denominaremos WINPRETOR ALFA. Por la gran complejidad de su código ha sido necesario destinar la mayor parte del tiempo a analizar el código fuente, por lo tanto se deja para futuros desarrollos la finalización de la aplicación. Pág. 8 Memoria 2. Contexto energético actual La humanidad requiere grandes cantidades de energía para poder cubrir sus necesidades más básicas, actualmente esas necesidades energéticas están cubiertas principalmente mediante el uso de combustibles fósiles. El uso de combustibles fósiles tienen dos problemas fundamentales; el primero, es el impacto medioambiental que produce su utilización; y el segundo es que por su naturaleza son ilimitados. Como hemos podido comprobar en anteriores ocasiones, la escasez de crudo produce consecuencias catastróficas en la economía mundial. Y esto es debido a que prácticamente todo lo que se produce en este planeta tiene un contravalor en oro negro. Por ejemplo, algo tan común como una vaca, tiene una equivalencia indirecta en barriles de crudo. Se ha demostrado que son necesarios tres barriles de crudo para poder desarrollar una vaca hasta su estado adulto, este computo tiene en cuenta la energía necesaria para procesar los piensos y la electricidad requerida para la iluminación, calefacción, maquinaria, etc. En la crisis del 74, la escasez de crudo fue producida por la coyuntura geopolítica del momento, y no precisamente por su escasez real, ya que las reservas de crudo eran baratas de explotar y había en abundancia. En el contexto actual, la escasez es consecuencia de un problema físico y muy real. Estamos llegando a los límites de capacidad de extracción de los pozos a nivel mundial. El crudo a corto plazo no se va a acabar, pero será cada vez más caro y difícil de extraer. El encarecimiento del crudo abre la puerta a nuevas oportunidades energéticas, como la fisión nuclear u o otras fuentes alternativas que hace unos años eran absolutamente inviables económicamente como por ejemplo la fusión nuclear. La fusión nuclear es la única fuente energética capaz de abastecer las necesidades mundiales en un futuro a largo plazo, ya que su rendimiento energético es el equivalente al de la fisión nuclear y no comporta sus problemas ecológicos. Pero, en contra tiene una gran complejidad tecnológica y requiere grandes inversiones en I+D antes de su posible explotación comercial. Muchos expertos empiezan a hablar de la economía del hidrógeno como principal candidato a substituir la actual economía del petroleo. En esta economía el medio de transporte de la energía dejaría de ser el petroleo para ser el hidrógeno. Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 9 El problema fundamental es como generar ese hidrógeno, puesto que en la naturaleza el hidrógeno como elemento independiente no existe ya que éste se oxida muy rápidamente y de una forma muy exotérmica produciendo uno de los compuestos más estables de la naturaleza: el agua. Para generar hidrógeno partimos del agua, iniciaríamos el proceso separando sus respectivos átomos de oxígeno e hidrógeno, para ello necesitamos gran cantidad de energía, precisamente la misma que nos cederá ese hidrógeno cuando lo volvamos a oxidar en el momento de su utilización final. Por lo tanto, para generar hidrógeno necesitaremos una fuente de energía muy potente y que genere la menor cantidad posible de residuos, esa fuente de energía es la fusión nuclear. El hidrógeno se generaría en centrales de fusión nuclear y se distribuiría para todo tipo de aplicaciones. El ciclo empieza con la extracción del Deuterio contenido en el agua de mar (un bien ilimitado) y acaba en la combustión de hidrógeno produciendo finalmente agua en forma de vapor, que volverá al mar a traves de los ciclos atmosféricos naturales. Un poder para la humanidad simplemente ilimitado. Pág. 16 Memoria Si éste es de suficiente magnitud, la partícula sale reflejada o simplemente frenada según la ecuación 4.4. Esta barrera magnética es ineficaz para partículas con µ pequeño o con velocidades muy altas, produciendo la pérdida sistemática de las partículas más energéticas que precisamente son las que más nos interesan. Las configuraciones abiertas tienen un problema, que cuando llegamos a temperaturas muy altas tienen tiempos de confinamiento de la energía muy bajos, impidiendo llegar a la ignición. No obstante, pueden ser útiles como generadores de intensos flujos de neutrones rápidos. Las configuraciones cerradas se basan en geometrías toroidales en las cuales las líneas de campo se cierran sobre si mismas. Mediante bobinas exteriores se genera un campo magnético toroidal que superpuesto con otro campo poloidal consigue confinar las partículas durante el tiempo suficiente para permitir la ignición. Dado que tenemos dos campos magnéticos superpuestos, las líneas de campo van enroscándose en el toro a medida que completan su trayectoria engendrando superficies anidadas que son conocidas como superficies magnéticas. En equilibrio, cada superficie magnética mantiene el mismo valor de densidad, presión y temperatura. F=-µgrad B (Ec. 4.4) Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 17 3.3. Reactores Tokamak En este trabajo nos centraremos en este tipo de reactores debido a que el software desarrollado se basa en éstos para la simulación. No obstante, existen otras líneas de investigación abiertas con prometedores resultados sobre otros reactores de configuración cerrada como los Stellator, que no trataremos en este trabajo. Tokamak son las siglas de toroidalnya kamera y magnetnaya katushka que significa cámara toroidal y bobina magnética. Como su nombre original en Ruso indica, la invención y el primer desarrollo de un Tokamak fue en la antigua URSS en los años 50. Actualmente, existen varios tokamaks de gran tamaño en el mundo como el TFTR en los EEUU, el JET de la UE y el Japonés JT60. A parte de estos existen multitud de reactores de tamaño pequeño y medio, orientados a la investigación en física de plasma y materiales. Como se muestra en la figura 4.4, el campo magnético helicoidal que confina el plasma es el resultado de la superposición de un campo magnético toriodal y otro poloidal. Campo magnético Toroidal: Se genera en un número finito de espiras que van rodeando la estructura generando el campo resultante que recorre todo el toro por su interior. Este campo por si solo produce un confinamiento intrínsecamente inestable dado que por la curvatura finita del toro, el campo magnético es más intenso en su parte más interior que en su parte más exterior. Como hemos visto anteriormente en la Ec. 4.5 un incremento en la intensidad del campo magnético genera una fuerza que en este caso repele las partículas hacia la parte más exterior del toro, estas partículas en su desplazamiento normal al campo magnético (drift) experimentan una fuerza de Lorentz perpendicular al campo y a su velocidad de desplazamiento, por lo tanto en el eje z (vertical). La fuerza de Lorentz depende del signo de cada partícula, por lo tanto las partículas positivas son acumuladas en la parte superior del toro y las negativas en la parte inferior. Esta situación acaba produciendo un campo eléctrico que llevará inevitablemente a la disrupción del confinamiento y al contacto del plasma con las paredes del reactor. La única solución a este grave problema pasa por generar otro campo magnético que haga ‘enroscarse’ las líneas de campo a medida que viajan por el toro para que así las partículas pasen el mismo tiempo en la parte más interior como en la más exterior evitando la formación del fatal campo eléctrico. Pág. 18 Memoria Fig. 4.4. Partes de un tokamak [4]. Campo magnético Poloidal: Este campo es necesario para estabilizar el plasma, el problema es como crearlo ya que como sabemos un campo magnético radial debe generarse por una corriente normal pasando por su interior. Así la única forma de generar este campo es mediante una corriente eléctrica en el seno del plasma, que recorre el toro como si de una espira se tratase. La forma clásica de generar esa corriente eléctrica en el plasma es mediante inducción electromagnética. A través de una bobina exterior al toro, utilizamos este como si fuera el secundario de un gran transformador (ver fig. 4.4). Este sistema es muy efectivo ya que además de generar la corriente eléctrica necesaria para estabilizar el plasma, conseguimos calentarlo ohmnicamente. Este calentamiento es efectivo sólo hasta temperaturas del orden de 1 Kev dado que a partir de aquí la resistividad del plasma es inferior a la del cobre y no conseguimos calentamientos ohmnicos significativos. Este sistema tiene el gran inconveniente de que para generar una corriente inducida en el plasma necesitamos un campo magnético variable en el tiempo (ya que es la variación del flujo magnético por una superficie cerrada por un conductor la que genera la corriente inducida en este), por lo tanto, la intensidad en el primario de nuestro transformador deberá crecer indefinidamente para que el reactor funcione de forma ininterrumpida. Como es evidente esto no es posible y tendremos que contentarnos con un reactor pulsante. Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 19 Se están estudiando formas de mantener corriente suficiente en el plasma sin tener que inducirla magnéticamente. Gracias a la corriente generada por la inyección de Neutros (NBI) y a la corriente de Bootstrap (generada en el plasma espontáneamente por los gradientes de densidad entre las diferentes superficies magnéticas y otras inestabilidades microscópicas) es posible mantener una corriente en este permitiendo su confinamiento. El problema radica en que estos plasmas no consiguen corrientes no inductivas lo suficientemente elevadas como para reemplazar la corriente inductiva en regímenes de potencia elevada. Pág. 20 Memoria 3.4. Balance de energía Un reactor funcionando en régimen estacionario debe cumplir la siguiente ecuación: Esta ecuación establece el equilibrio térmico global del plasma igualando las fuentes de energía con sus pérdidas. En la parte izquierda de la ecuación nos encontramos las fuentes de energía y en la derecha las pérdidas. Fuentes: Pα: Es la potencia depositada en el plasma por las partículas Alfa, resultado de las reacciones de fusión. Esta energía es el 20% del total producido por las reacciones y el 80% restante desaparece del plasma en forma de neutrones a 14 MeV y por ello no aparece en el balance. POH: Es la potencia generada sobre el plasma mediante la inducción de la corriente Toroidal. La POH depende, de una forma muy importante, de la temperatura del plasma ya que, la resistividad de éste baja con su temperatura. A partir de 1 KeV su contribución al calentamiento del plasma es mínima. Padd: Son las diferentes fuentes de energía adicionales tales como, haces de neutros o radiaciones electromagnéticas resonantes. Pérdidas: PB son las pérdidas producidas por la radiación de Bremsstrahlung que tiene su origen en los electrones cuando cambian de dirección en presencia de núcleos atómicos. Esta radiación aumenta con la temperatura de los electrones a razón de Te1/2 y de una forma crucial con Zeff. Es importantísimo minimizar la cantidad de impurezas en el plasma y las que hayan, que sean de una Z lo más pequeña posible. Sus emisiones están en la frecuencia de los rayos X. (Ec. 4.5) Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 21 PSyn son las pérdidas del efecto Sincrotón, por el giro constante de los electrones en sus trayectorias helicoidales. Como sabemos, cualquier partícula cargada cuando cambia de trayectoria emite radiación en la frecuencia correspondiente al cambio de energía cinética sufrido. La radiación producida por el efecto Sincrotón tiene una longitud de onda alrededor del milímetro, siendo para estas longitudes de onda el plasma muy opaco y las paredes del reactor muy reflectantes. Una gran parte de esta energía es reabsorbida por el mismo plasma. Hay que decir que empieza a cobrar importancia a temperaturas muy altas. (Ec. 4.6) (Ec. 4.7) Pág. 22 Memoria Pcon: Es la cantidad de energía transferida desde el centro a la frontera exterior del plasma , por procesos de conducción y convección. Su simulación es una de las más complejas porque depende de los fenómenos de transporte radial de energía a través de las diferentes superficies magnéticas. Un forma simple de modelar este coeficiente es mediante el tiempo de confinamiento de la energía (Ec.4.8). La potencia perdida por convección es resultado de la división de la energía total contenida en el plasma por el tiempo de confinamiento. (Ec. 4.8) Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 23 Fig. 5.1. Estructura de la simulación PRETOR. 4. Modelo PRETOR Como he comentado anteriormente, para la realización de este trabajo me he basado en el código fuente del simulador PRETOR PC que incorpora el código de modelización de plasmas PRETOR. Este código permite simular radialmente las variaciones de las diferentes magnitudes físicas del plasma en función del tiempo. La simulación es unidimensional y discretizada en diferentes nodos que vendrían a representar las superficies magnéticas. El código está dividido en diferentes módulos, el módulo principal es el encargado de resolver los sistemas de ecuaciones no-lineales acopladas, correspondientes al transporte de la energía y partículas. Estos cálculos suelen plasmarse en coeficientes dinámicos de transporte que son calculados mediante el modelo Heurístico Rebut-Llalia-Watkins [2]. El modelo PRETOR también permite añadir factores adicionales a la simulación, como formas auxiliares de calentamiento por haces de neutros (NBI) y por radiofrecuencia (ECRH, ICRH). Nos permite tener el inventario de partículas en cada instante de la simulación, incluyendo la alimentación de combustible, la producción de cenizas y la densidad de impurezas. Parámetros Iniciales del plasma Inicialización variables de entorno Parámetros de plasma a t Parámetros de plasma a t+n t=0 Código PRETOR t=t+n Pág. 24 Memoria Existe la posibilidad de simular diferentes efectos MHD como los dientes de sierra o el balloning. Como podemos ver en la figura 5.1, el código PRETOR está compuesto por un conjunto de algoritmos que a partir de un estado inicial de las variables del plasma, van evolucionando en el tiempo. Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 25 4.1. Estructura del Código PRETOR PC Pretor PC es una programa informático creado por el Dr. Boucher en el año 1992. Dado que por ese entonces Windows era una rareza, está diseñado en entorno MSDOS y en lenguaje Basic. Al inicio del proyecto dediqué innumerables horas a intentar compilar los archivos Basic, probé todos los compiladores disponibles en el mercado. Todo esfuerzo fue inútil, el código fuente está lleno de errores de sintaxis y una vez corregidos, recibía errores por parte del compilador del tipo ‘Out of Memory’. El resultado ha sido que no he podido compilar ni ejecutar paso por paso el programa, ello ha aumentado muchísimo la comprensión de código fuente que para rematar no incluía ni un solo comentario. 4.1.1. Ficheros El código fuente consta de 8 archivos .BAS y 4 archivos .BI. En los archivos BAS se encuentra el código de las rutinas y en los archivos BI declaraciones de variables globales y estructuras de datos. Pág. 32 Memoria Nombre DESCRIPCIÓN COORD Es una estructura que permite definir puntos en coordenadas toroidales, es utilizado por la matriz PointRZ(0, NbB) para definir la trayectoria de los pellets y de la inyección de neutros. INFO Se aplica al vector global inf(240) que contiene información de todas las gráficas disponibles en memoria y en pantalla. La variable Kind es muy importante ya que especifica el tipo de información que contienen el resto de variables. 0= Radial profile. 1= Exp. radial profile. 2= Text. 3=Time traces at given radii. 4=Time traces. 5=Empty. Otra variable de vital importancia es Curb que indica el índice de la matriz PHYS al que corresponde la variable de la gráfica. Echmin/Echmax corresponden al Máximo y mínimo en gráfica. El resto de parámetros tendrán información dependiendo del valor de Kind: KIND: kind name curb Decc Emin Emax Radius ETime -2 X -1 X X 0 X X X X X X 2 X X 3 X X X X X X 4 X X X X X 5 X Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 33 Fig. 5.6. Descripción de cada estructura de datos en el código fuente PRETOR PC. INFTIMEXP Se aplica al vector global InfEx(16) que contiene información necesaria para procesar las trazas de tiempo experimentales como el inicio temporal el final el máximo y el mínimo. EVREC Estructura que se aplica al vector global EventR(MaxEvent) que contiene información sobre los eventos temporales. Es de suma importancia para la simulación. DEPOSITION Estructura utilizada en multitud de partes del código para representar el valor de una variable a lo largo de una trayectoria, por ejemplo es utilizada para definir el grado de transferencia de energía cinética de los haces de neutros a medida que van penetrando en el plasma. ANIEXP Estructura utilizada a través de la variable de ámbito global anime para almacenar la información de cabecera de todas las trazas de tiempo experimentales cargadas. Pág. 34 Memoria 4.1.5. Variables Globales significativas PRETOR guarda el estado de la simulación en multitud de variables globales, que son consultadas y modificadas por las funciones de simulación, interficie de usuario y de ficheros. La mayor parte de la información de la simulación está contenida en la matriz de variables de estado del plasma denominada PHYS(n, m). Esta matriz es de n filas y m columnas, siendo n el número de nodos que modelan el plasma y m cada una de las variables. A continuación incluyo una tabla con todas las posiciones de la matriz y la variable a que hacen referencia: Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 35 Fig. 5.7. Índice de variables de la matriz PHYS. Te = 1 j boot = 16 DpartI = 31 NeutrIC = 46 VOL = 61 Ti = 2 grdCrit = 17 VpinchI = 32 NeutrIH = 47 VD = 62 ne = 3 KhiE = 18 DpartImp = 33 NeutrP = 48 ni = 4 KhiI = 19 vpinchImp = 34Emiss = 49 np = 5 DpartE = 20 Baloo = 35 FLUX = 50 nion = 6 SELEC = 21 grdP = 36 Z2 = 51 zf = 7 SIONS = 22 KhiBal = 37 WE = 52 z = 8 PBrems = 23 DPartZ = 38 WI = 53 j = 9 RpadHF = 24 vpinchZ = 39 meqPast = 54 q = 10 RNBI = 25 SnbPart = 40 TeBS = 55 e = 11 cyclo = 26 RpadnbE = 41 TiBS = 56 Bt = 12 PalfE = 27 RpadnbI = 42 PhiE = 57 Bpol = 13 PalfI = 28 SHEPart = 43 Test = 58 It = 14 meq = 29 Empty = 44 PRFelec = 59 eta = 15 poh = 30 NeutrI = 45 PRFion = 60 Pág. 36 Memoria Fig. 5.8. Ejemplo de matriz del contenido de la matriz PHYS. En la figura 5.8 podemos ver un ejemplo de esta matriz en un instante de la simulación. Id N1 N2 N3 N4 ... N(MAX) 1 Te 6,29E-01 9,75E-01 1,62E+00 2,39E+00 ... 8,28E+00 2 Ti 1,19E+00 1,19E+00 1,75E+00 2,48E+00 ... 7,87E+00 3 Ne 4,13E+00 5,76E+00 6,25E+00 6,61E+00 ... 7,73E+00 4 Ni 4,10E+00 5,68E+00 6,16E+00 6,49E+00 ... 7,51E+00 5 Np 1,34E-02 4,12E-02 4,79E-02 5,57E-02 ... 1,10E-01 ... ... ... ... ... ... ... ... 58 Test 1,01E+00 1,01E+00 1,02E+00 1,02E+00 ... 1,03E+00 Otra vector de vital importancia es PVAL(1-160) que contiene todos los parámetros recuperados del fichero .PAR desde la función LoadDat. Estos valores son utilizados constantemente para realizar la simulación y son los que podemos modificar en las diferentes pantallas de PRETOR PC. El nombre de cada uno de los parámetros está contenido en el vector Ptitle (1-160). Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 37 Fig. 5.9. Ejemplo de contenido de los vectores Pval y Ptitle. Indice Nombre V alo r 2 Speed 0,005 3 Display 1 4 Mesh 14 11 -Curb scaling: 0 12 Pannel Number 10 13 Curb number 1 14 Lower limit 0 15 Upper limit 20 16 Radius (m) 0 ... ... ... 160 P CD / j : 0 A continuación incluyo un resumen de las variables globales que he considerado más significativas: Nombre de variable Dimesiones Subsistema Descripcion Machine 0 Simulación Identifica el nombre de la máquina que se va a utilizar. Cada una de las máquinas están modelizadas siguiendo unos parámetros diferentes en el código que no se pueden modificar. La máquina por defecto es JET pero puede ser TFTR, ITER, RTP, DIIID, TORE SUPRA, START. INF(240) 1 Gráfico Vector de la estructura INFO donde se almacena la información sobre todas las variables gráficas. Son destacables las variables: Curb que indica la posición en la matriz PHYS de la variable física y Iname que indica el nombre de la gráfica. Pág. 38 Memoria Fig. 5.10. Algunas de las variables más importantes de Pretor PC con una descripción. PanTle(1-24) 1 Gráfico Nombre de las 8 pantallas de PRETOR PC Ptitle(1-160) 1 Gráfico Nombre de todas las variables del reactor. Puede llegar hasta 240 pero el programa sólo utiliza 160. Está ligada a la variables donde están almacenados los valores de cada variable PVAL(1-160) PVAL(1-160) 1 Simulación Valor de las variables globales del reactor. Se inicializan con un valor al inicio de la aplicación en Winpretor en la función CargarPanels. Luego son cargadas desde el fichero .PAR correspondiente al escenario seleccionado en la función LoadDat. PHYS(62, ADim) 2 Simulación Es la matriz más importante de todo el programa porque es donde se almacenan todas las variables de estado del plasma en cada uno de sus nodos. La primera dimensión de la matriz, corresponde al código numérico asociado a cada una de las variables, la segunda dimensión tiene una longitud impuesta por la variable global ‘max’ que viene dada en la primera fila del fichero .PAR y significa el numero de superficies magnéticas. Physname(1-58) 1 Gráfico Nombre de cada una de las variables contenidas en la matriz PHYS. MAX 0 Simulación Es el número de nodos que modelan cada una de las superficies magnéticas en el perfil radial del plasma. Tax 0 Simulación Tiempo en segundos de la simulación. dT 0 Simulación Diferencial de tiempo entre dos pasos de la simulación, utilizado en todas la ecuaciones diferenciales respecto al tiempo. Es variable ya que se ajusta dinámicamente a las condiciones del plasma. Simplificando mucho podríamos decir que es inversamente proporcional a las variaciones de energía del plasma. De esta forma afina más en los transitorios rápidos. Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 39 Fig. 5.11. Ciclo de simulación simplificado. 4.1.6. Diagrama de flujo simplificado del motor de simulación Una vez todas las matrices y variables de plasma han sido cargadas en memoria mediante las diversas funciones de carga de ficheros (principalmente LoadDat), el programa entra en el bucle de simulación que está contenido en la función Main desde la cual para cada paso de la simulación se accede a la función Keyboard, que es la que controla la actividad del teclado. En caso de producirse una pulsación de teclado ésta es detectada en la función Keyboard que actuará de la forma correspondiente. A continuación explicaré los pasos más importantes en un ciclo de simulación a medida que se van produciendo, evidentemente sólo explico lo más importante dado que la complejidad de la aplicación es muy grande. Main (Bucle de simulación) 1. Recalculo dt a aplicar 2. Proceso de eventos 3. Rutina MainPhys 4. Calculo del fenómeno Sawtooth 5. Keyboard MainPhys 1. Verifico Pellets 2. Cargo variables temporales 3. Comprueba calentamiento RF 4. Comprueba rampas de todo tipo 5. Comprueba haz de neutros 6. Calculo de perdidas por radiación Cyclotrónica 7. Calculo potencia Alfa 8. Calculo por Bremsstrahlung, la densidad de potencia de los electrones e iones, 9. la potencia ohmnica y la variable Test 10. Calculo coeficientes de difusión térmicos de electrones e iones 11. Calculo coeficientes de difusión de partículas de electrones e iones 12. Recalculo delinventario de partículas 13. Recalculo temp. De los iones y electrones 14. Calculo de la energia contenida en el plasma 15. Corriente de BootStap 16. Calculo campo magnético poloidal 17. Calculo densidad de corriente, campo magnético toroidal, factor de seguridad. 18. Fenómeno Balooning 19. Calculo del flujo total de energía entre superficies magnéticas 20. Comprobación del fenómeno Sawtooth Pág. 40 Memoria • Recalcula el tiempo para todas las gráficas de experimentos. • Procesa los eventos asignados a ese instante de tiempo • Recalcula el dt de tiempo a aplicar en el siguiente paso de la simulación. El cálculo se realiza teniendo en cuenta el incremento de energía en el plasma en el paso interior, si este incremento ha sido muy grande, el programa entiende que estamos ante un transitorio brusco y disminuye el dt para que la simulación sea más precisa. Por lo tanto, dt tiene que estar en unos límites (dtmax, dtmin) pero puede variar dependiendo de las condiciones dinámicas del plasma. • Entra en la función MainPhys que concentra todas las operaciones físicas del plasma. • Verifica que no haya entrado un Pellet mediante un evento, esto se comprueba mediante el valor de Npellet, si este es mayor que 0 entra en la rutina de proceso de los pellets llamada Pellet. • Carga diferentes filas de la matriz de estado del plasma en matrices temporales para poder trabajar con ellas de una forma más sencilla y rápida. • Comprueba el calentamiento por radiofrecuencia (PVAL(160) <> 0), si está activado actualiza las variables necesarias para que se tenga en cuenta en el recálculo de los balances energéticos mediante la rutina ProfHF. • Comprueba que no esté ninguna rampa activada, en caso afirmativo entra en la función RAMP que es la encargada de actualizar el valor sobre el que actúa la rampa y calcular su avance hasta el valor máximo en cada uno de los siguientes ciclos del bucle. • Verifica que no esté activado el haz de neutros, en caso afirmativo entra en la rutina de cálculo de neutros NBIProfile. • Calcula las pérdidas por radiación ciclotrónica en el plasma mediante la función Pcyclo. Almacena el resultado en el vector CPcyc() de ámbito global que nos indica la potencia emitida por cada uno de los nodos. • Calcula la potencia alfa para cada uno de los nodos mediante una función llamada ALFA. Si la ni (densidad iónica del combustible) es mayor que cero y la Ti (temperatura electrónica) es mayor que 1 Kev(*) se considera que se está Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 41 produciendo una cantidad significativa de potencia alfa y se calcula para cada uno de los nodos. * Recordemos que a 1 KeV la potencia Alfa es unas 300 veces menor que a 13,5 KeV. Por lo tanto por debajo de esta temperatura es correcta la suposición de igualarla a cero. • En la función HeatSources calcula para cada una de las superficies magnéticas o nodos: el Flujo de electrones, las pérdidas de radiación por Bremsstrahlung, la densidad de potencia de los electrones e iones, la potencia ohmnica y la variable Test utilizada por el programa para decidir el momento de transición* a modo H de confinamiento magnético mejorado. * Si el valor de Test en la 2ª superficie magnética exterior es mayor que PVAL(144) el programa decide que estamos en modo H. • En la función CoefHeat se calculan por cada superficie magnética los coeficientes de difusión térmica de los iones y los electrones en t+1. También se calcula el parámetro GrdCrit, su valor determina el momento donde empieza el transporte de la energía anómalo, mientras tanto se aplica el modelo neoclásico. • La función CoefPart calcula por cada superficie magnética los coeficientes de difusión de partículas en t+1, de impurezas y sus respectivos coeficientes de convección. • Recalcula el inventario de partículas t+1 en el plasma por superficie magnética incluyendo la generación de He, las impurezas, las posibles inyecciones de combustible en forma de pellets, las pérdidas por combustión y escape del confinamiento. En este proceso interviene la rutina PartFlow ya que recalcula todos los flujos de las diferentes partículas entre superficies magnéticas. • Recalcula las temperaturas de iones y electrones en el instante t+1 y en cada superficie magnética a través de la rutina HeatDiff que se aplica a iones y electrones indistintamente. • Calcula la energía almacenada por los electrones e iones del plasma y en cada superficie magnética en el instante t+1 con la rutina Energ que se aplica a iones y electrones indistintamente. • Se Calcula la corriente de BootStrap o corriente autogenerada por el plasma a t+1 en cada superficie magnética mediante la rutina CalcEta • Cálculo del campo magnético Poloidal a t+1 mediante la rutina FieldDiff en cada superficie magnética. Pág. 48 Memoria Las pantallas secundarias contienen la siguiente información: Perfiles de plasma de hasta 15 variables simultáneas en una misma gráfica comparativa. Evoluciones temporales de variables del plasma pudiendo contener un nodo, la media de todos o simplemente todos a la vez. Variables integradas globales de la simulación como potencias totales. Balances energéticos. Perfiles geométricos del plasma con sus características físicas principales. Visualización de parámetros del fichero .PAR. Visualización de todas las variables de la matriz PHYS(n,m). Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 49 Fig. 6.3. Botón de abrir proyecto en la pantalla principal de WinPretor Alfa. 5.2.2. Módulo de carga de ficheros Sólo permite la carga de ficheros .PAR y .DAT. Se ha desactivado la posibilidad de carga de otros tipos de ficheros como por ejemplo las trazas de tiempo. Evidentemente, se ha mantenido la compatibilidad de los ficheros con la versión DOS y LINUX. Además carga un nuevo tipo de fichero llamado .WPR que es el encargado de guardar y recuperar el estado de la pantalla. Es decir, almacena la información relativa a la posición de cada ventana y su actividad. A continuación veremos paso a paso como se abre un proyecto desde WinPretor Alfa. Botón de carga de proyecto. Pág. 50 Memoria Fig. 6.4. Pantalla de ficheros de WinPretor Alfa. Desde la pantalla (Fig 6.4) primero tenemos que buscar el proyecto que queremos cargar, una vez seleccionado podemos ver la máquina que utiliza, el nombre del escenario y el número de nodos que utiliza. En la parte inferior de esta pantalla podemos especificar las opciones de interpolación nodal, esta es una nueva funcionalidad que nos permite multiplicar el número de nodos por 2 ó por 4. El aplicativo parte de los nodos iniciales y aplica una interpolación lineal para encontrar los nuevos nodos, a partir de este momento la simulación se desarrolla con una mayor precisión. Buscar proyecto a abrir. Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 51 Fig. 6.5. Resultado de las diferentes interpolaciones internodales. 14 Nodos 28 Nodos 56 Nodos 2X 4X El siguiente paso es hacer clic en el botón de Abrir, en este momento el programa accede a los ficheros .PAR, .DAT, WPR , los carga en memoria y los procesa. En la pantalla principal aparece el proyecto con la última configuración de pantalla. Pág. 52 Memoria Fig. 6.6. Proyecto cargado en el entorno MDI. Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 53 Fig. 6.7. Opciones de simulación (paso por paso o a saltos de 10 seg). Fig. 6.8. Tiempo total de la simulación y el dt aplicado. 5.2.3. Módulo de simulación simplificado En esta versión de WinPretor sólo se ha incluido el módulo de simulación simplificado que permite la evolución de un estado de plasma inicial a partir de los perfiles radiales cargados desde el fichero .DAT y los parámetros generales de la simulación .PAR. A partir de este momento no se puede modificar nada y la simulación avanza hasta que la paremos. La simulación puede visualizarse paso por paso o en intervalos de 10 seg. En la parte inferior de la pantalla se muestra continuamente el tiempo total de la simulación y el dt aplicado en el último paso, este dt es variable en el tiempo. En la parte superior izquierda está la información del escenario cargado en memoria, como pueden ser el número de nodos, la máquina utilizada y el nombre del escenario. Pág. 54 Memoria Fig. 6.9. Panel de información general de la simulación. Fig. 6.10. Pantalla de visualización de los parámetros de la simulación. 5.2.4. Pantalla de visualización de parámetros Los parámetros son las variables fijas del reactor y de la simulación, que se encuentran en el fichero .PAR, estas pueden ser consultadas desde la pantalla de visualización de variables del reactor. En la versión actual sólo se pueden visualizar y no se pueden modificar. Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 55 Fig. 6.11. Pantalla de visualización de variables del plasma. 5.2.5. Módulo de visualización de variables del plasma Es el módulo más importante ya que de este parten el módulo de visualización de perfiles y evoluciones temporales de variables. Básicamente, es una pantalla en la que se muestran, en forma de columnas, los nodos del plasma, y en forma de filas, todas las variables de estado. De esta forma, en una misma pantalla, tenemos toda la información del plasma en un momento dado. Seleccionando una fila, podemos ver en la parte inferior de la pantalla una representación gráfica de la variable relacionada. Es posible seleccionar varias filas y ver diferentes gráficas simultáneamente hasta un máximo de 15. Desde esta pantalla se puede exportar la matriz de parámetros a Excel desde el botón indicado en la figura 6.11. Pág. 56 Memoria Fig. 6.12. Pantalla de perfil radial de plasma independiente mostrando visualización gráfica. Haciendo clic con el botón derecho del ratón sobre una fila es posible acceder al menú de acciones, que nos permite generar una ventana de perfiles radiales independiente o una ventana de evolución temporal. De esta forma, a partir de la pantalla de variables del plasma, vamos creando todas las visualizaciones que queramos sobre las variables que más nos interese controlar. En los dos siguientes apartados veremos como funcionan las ventanas de evoluciones temporales y perfiles de plasma. 5.2.6. Módulo de perfiles de plasma Este módulo permite generar los perfiles del plasma que más nos interesen a partir del módulo de variables del plasma. Las pantallas independientes de perfiles irán evolucionando con la simulación a medida que esta avance. Selección del tipo gráfica. Selección del modo de visualización, gráfica o datos. Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 57 Fig. 6.13. Diferentes tipos de gráficas. Como podemos ver en la Fig. 6.12, desde cada venta de perfiles podemos seleccionar el tipo de gráfica y cambiar la visualización de gráfica a datos (Fig. 6.13). Pág. 64 Memoria Conclusiones Con este proyecto se ha desarrollado la versión WINPRETOR ALFA en la cual se han implementado la mayor parte de las rutinas de cálculo de PRETOR y se ha desarrollado un nuevo entorno gráfico agradable y fácil de utilizar. Todos los objetivos fijados al inicio del proyecto han sido superados ampliamente, ya que en un inicio pensábamos que como mucho llegaríamos a documentar el código PRETOR y finalmente hemos acabado desarrollando una versión funcional. Este proyecto pretende ser el inicio de uno más ambicioso, que comprendería el desarrollo, en su totalidad, de un simulador de fusión nuclear que incorporaría los últimos algoritmos de transporte de plasma y utilizaría la más moderna tecnología de computación distribuida para poder realizar simulaciones con una precisión sin precedentes. Personalmente ha sido muy satisfactorio el desarrollo de este Proyecto y es mi intención darle continuidad durante el doctorado. Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 65 Agradecimientos En primer lugar me gustaría agradecer la ayuda y conocimientos prestados por mi director de proyecto, el Dr. Guillem Cortes. También me gustaría agradecer la inestimable ayuda de Susana en la redacción y la ortografía. Y por supuesto agradecer a mi familia el soporte incondicional y la paciencia prestada durante todos estos años de carrera. Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 67 Bibliografía Referencias bibliográficas [1] SECCIÓN DE FÍSICA E INGENIERÍA NUCLEAR (UPC) Javier Dies, Ferran Albajar, Joan Fontanet, “Utilización de un simulador de un reactor de fusión nuclear tipo Tokamak para fines docentes”. Barcelona, Octubre 1998 [2] DOMINIQUE BOUCHER “Etude et modélisation du transport de l’ energie et des particules dans un plasma de fusion thermonucléaire contrôlée. Application au transport anormal et aux conditions de fonctionnement du Tokamak.” París, Septiembre 1992. [3] J. Raeder et al. “Controlled nuclear fusion. Fundamentals of its utilization for energy supply” Wiley & Son, New Cork, 1986. [4] De Barbieri, O., “Sincrotron Emissions from High Temperatura Plasmas” Report EURCEA-FC-1035, 1980. [5] FOM-Institute for Plasma Physics Association EURATOM-FOM. “The Fusion Road Show Presentation” 2003 [6] San Diego Research Institute, “Fusión, Creating a Star on Herat”, 2001 Pág. 68 Memoria Bibliografía complementaria: [1] Soft Conference, “How far is a Fusion Reactor from an Experimental Reactor?” Madrid, 2000. [2] 28th EPS Conference on Contr. Fusion and Plasma Phys. Funchal, “SAWTOOTH STABILISATION BY NEUTRAL BEAM-INJECTED FAST IONS IN JET” Junio 2001 [3] IAEA “Summary of the Iter Final Design Report” Viena 2001 [4] K. Bockasten, R. Hallin, S.I. Herlitz, Viena International Atomic Energy Agency, “Controlled thermonuclear fusion, research” Viena 1961 [5] Duchs, D. F., “Recent steps towards a controlled thermonuclear fusion reactor with results from the JET tokamak device / D. F.”, 1992 [6] Página Web: http://www.iter.org Optimización del código del simulador de reactores de fusión nuclear PRETOR Pág. 69