html5-img
1 / 111

Sistemas de Ecuaciones Lineales

Sistemas de Ecuaciones Lineales. Glen D. Rodríguez R. Algoritmos Paralelos Basado en material de J. Demmel. Eliminación Gaussiana y matrices densas. Algebra lineal densa “Gaussian Elimination” (GE) para resolver Ax=b Optimizar GE para máquinas secuenciales con caché

halona
Télécharger la présentation

Sistemas de Ecuaciones Lineales

An Image/Link below is provided (as is) to download presentation Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript

  1. Sistemas de Ecuaciones Lineales Glen D. Rodríguez R. Algoritmos Paralelos Basado en material de J. Demmel

  2. Eliminación Gaussiana y matrices densas • Algebra lineal densa • “Gaussian Elimination” (GE) para resolver Ax=b • Optimizar GE para máquinas secuenciales con caché • Usando multiplicación matriz x matriz (BLAS) • LAPACK • Partición de la Data en comp.paralelas • GE paralelo • ScaLAPACK • Problemas de Eigenvalue

  3. Sca/LAPACK Overview

  4. Exitos usando Sca/LAPACK • Muy usado: • Adoptado por Mathworks, Cray, Fujitsu, HP, IBM, IMSL, NAG, NEC, SGI, … • >56M visitas en la web de Netlib (incl. CLAPACK, LAPACK95) • Descubrimientos hechos a través de la solución de sistemas de ec. densos • Artículo en Nature sobre el universo plano usó ScaLAPACK • Otros en Physics Review B también lo usaron • 1998 Premio Gordon Bell • www.nersc.gov/news/reports/newNERSCresults050703.pdf Cosmic Microwave Background Analysis, BOOMERanG collaboration, MADCAP code (Apr. 27, 2000). ScaLAPACK

  5. Usos de Solución de sistemas de ecuaciones lineales • Resistencia de estructuras (ing. civil) • Stress, tensiones en piezas mecánicas • Electromagnetismo, electrodinámica • Circuitos eléctricos • Ing. de tráfico. • DSP, óptica (de-convolución) • Economía y finanzas.

  6. Motivación (1) 3 problemas básicos de Alg. Lineal • Ecuaciones Lineales: Resolver x en Ax=b • Mínimos cuadrados: Hallar x que minimiza ||r||2S ri2donde r=Ax-b • Estadística: Ajustar datos a funciones simples 3a. Eigenvalues: Hallar l, x tal que Ax = l x • Análisis de vibraciones (terremotos, circuitos) 3b. Singular Value Decomposition: ATAx=2x • Ajuste de data, Obtención de información Muchas variaciones para diferentes estructuras de A • Simétrica, definida positiva, en bandas, …

  7. Motivación (2) • Por qué A densa, en vez de A dispersa? • Muchas matrices grandes son dispersas, pero … • Alg. Densos son más fáciles de entender • Algunas aplicaciones producen matrices densas grandes • LINPACK Benchmark (www.top500.org) • “Qué tan rápida es tu computadora?” = “Qué tan rápido resuelves Ax=b denso?” • Alg. de matriz dispersa grande con frecuencia producen menores (pero aún grandes) problemas densos.

  8. Records en la solución de Sist. Densos (2006) www.netlib.org Gigaflops Máquina n=100 n=1000 Any n Pico IBM BlueGene/L 281K 367K (131K procs) (281 Teraflops) (n=1.8M) NEC SX 8 (8 proc, 2 GHz) 75.1 128 (1 proc, 2 GHz) 2.2 15.0 16 … Palm Pilot III .00000169 (1.69 Kiloflops)

  9. Records en la solución de Sist. Densos (2011/2) www.netlib.org Performance (Gflops) , Problem size on 4 nodes. 4 AMD Athlon K7 500 Mhz (256 Mb) – (2x) 100 Mbs Switched – 2 NICs per node (channel bonding) 8 Duals Intel PIII 550 Mhz (512 Mb) - Myrinet Compaq 64 nodes (4 ev67 667 Mhz processors / node) AlphaServer SC Performance (Gflops), Problem size on 8- and 16-processors grids.

  10. Eliminación Gaussiana (GE) en Ax=b • Añadir múltiplos de cada fila a las demás para volver A en triangular sup. • Resolver el sist.triang. resultante Ux = c por sustitución … para cada columna i … “cerear” debajo de la diagonal añadiendo múltiplos de la fila i a las de abajo for i = 1 to n-1 … para cada fila j debajo de la fila i for j = i+1 to n … añadir un multiplo de la fila i a la fila j tmp = A(j,i); for k = i to n A(j,k) = A(j,k) - (tmp/A(i,i)) * A(i,k) … 0 . . . 0 0 . . . 0 0 . . . 0 0 . . . 0 0 . . . 0 0 . . . 0 0 . . . 0 0 . 0 0 . 0 0 0 0 After i=1 After i=2 After i=3 After i=n-1

  11. Refinando el algoritmo GE (1) … para cada columna i … “cerear” debajo de la diagonal añadiendo múltiplos de fila i a las de abajo for i = 1 to n-1 … para cada fila j debajo de fila i for j = i+1 to n … añadir un múltiplo de la fila i a la fila j tmp = A(j,i); for k = i to n A(j,k) = A(j,k) - (tmp/A(i,i)) * A(i,k) • Versión inicial • Elimina computación de constante tmp/A(i,i) de loop interno for i = 1 to n-1 for j = i+1 to n m = A(j,i)/A(i,i) for k = i to n A(j,k) = A(j,k) - m * A(i,k) m

  12. Refinando el algoritmo GE (2) for i = 1 to n-1 for j = i+1 to n m = A(j,i)/A(i,i) for k = i to n A(j,k) = A(j,k) - m * A(i,k) • Última versión • No computa lo que ya se conoce: ceros bajo la diagonal en la columna i for i = 1 to n-1 for j = i+1 to n m = A(j,i)/A(i,i) for k = i+1 to n A(j,k) = A(j,k) - m * A(i,k) m No calcula los ceros

  13. Refinando el algoritmo GE (3) for i = 1 to n-1 for j = i+1 to n m = A(j,i)/A(i,i) for k = i+1 to n A(j,k) = A(j,k) - m * A(i,k) • Última versión • Guarda los factores m bajo la diagonal en el espacio donde irían los ceros para uso posterior for i = 1 to n-1 for j = i+1 to n A(j,i) = A(j,i)/A(i,i) for k = i+1 to n A(j,k) = A(j,k) - A(j,i) * A(i,k) m Guarda m aquí

  14. Refinando el algoritmo GE (4) for i = 1 to n-1 for j = i+1 to n A(j,i) = A(j,i)/A(i,i) for k = i+1 to n A(j,k) = A(j,k) - A(j,i) * A(i,k) • Última versión • Partir el loop for i = 1 to n-1 for j = i+1 to n A(j,i) = A(j,i)/A(i,i) for j = i+1 to n for k = i+1 to n A(j,k) = A(j,k) - A(j,i) * A(i,k) Guardar todos los m aquí antes de modificar el resto de la matriz

  15. Refinando el algoritmo GE (5) for i = 1 to n-1 for j = i+1 to n A(j,i) = A(j,i)/A(i,i) for j = i+1 to n for k = i+1 to n A(j,k) = A(j,k) - A(j,i) * A(i,k) • Última versión • Rápido usando operaciones de matrices (BLAS) for i = 1 to n-1 A(i+1:n,i) = A(i+1:n,i) * ( 1 / A(i,i) ) A(i+1:n,i+1:n) = A(i+1:n , i+1:n ) - A(i+1:n , i) * A(i , i+1:n)

  16. Qué computa en realidad la GE? for i = 1 to n-1 A(i+1:n,i) = A(i+1:n,i) / A(i,i) A(i+1:n,i+1:n) = A(i+1:n , i+1:n ) - A(i+1:n , i) * A(i , i+1:n) • Llamemos a la matriz triangular inferior de multiplicadores m como M, y hagamos L = I+M • Llamemos al triangulo superior de la matriz final como U • Lema (Factorización LU): Si el algoritmo anterior termina (no hay división por cero) entonces A = L*U • Solución de A*x=b usando GE • Factorizar A = L*U usando GE (costo = 2/3 n3 flops) • Resolver L*y = b para y, usando substitución (costo = n2 flops) • Resolver U*x = y para x, usando substitución (costo = n2 flops) • Por lo tanto A*x = (L*U)*x = L*(U*x) = L*y = b como se desea

  17. Problemas con el alg. GE básico • Qué pasa si algún A(i,i) es cero? O es muy pequeño? • Resultado inválido (div.por cero), o “inestable”,  se necesita pivotear • Se computan operaciones BLAS 1 o BLAS 2, pero sabemos que las BLAS 3 (mult. de matrices) es más rápida (clase anterior…) for i = 1 to n-1 A(i+1:n,i) = A(i+1:n,i) / A(i,i) … BLAS 1 (escala un vector) A(i+1:n,i+1:n) = A(i+1:n , i+1:n ) … BLAS 2 (update en rango 1) - A(i+1:n , i) * A(i , i+1:n) Pico BLAS 3 BLAS 2 BLAS 1

  18. Pivotear en GE • A = [ 0 1 ] falla por que no se puede dividir entre A(1,1)=0 • [ 1 0 ] • Pero resolver Ax=b es facilísimo! • Cuando elemento diagonal A(i,i) es diminuto (no sólo cero), el algoritmo puede terminar pero dar una respuesta completamente errada • Inestabilidad numérica • La causa es el error de redondeo a punto flotante • Cura: usar Pivot (intercambiar filas de A) tal que A(i,i) sea grande

  19. GE con Pivot Parcial (GEPP) • Pivoteo Parcial: intercambiar filas tal que A(i,i) es la mayor en su columna for i = 1 to n-1 Hallar y guardar k donde |A(k,i)| = max{i <= j <= n} |A(j,i)| … o sea, el mayor element en el resto de la columna i if |A(k,i)| = 0 exit con advertencia que A es singular, o casi singular elseif k != i intercambiar filas i , k de A end if A(i+1:n,i) = A(i+1:n,i) / A(i,i) … cada cociente cae en [-1,1] A(i+1:n,i+1:n) = A(i+1:n , i+1:n ) - A(i+1:n , i) * A(i , i+1:n) • Lema: este algoritmo computa A = P*L*U, donde P es una matriz de permutación. • Es numéricamente estable en la práctica • Para más detalles ver código LAPACK en • http://www.netlib.org/lapack/single/sgetf2.f

  20. Problemas con el alg. GE básico • Qué pasa si A(i,i) es cero? O es muy pequeño? • Resultado inválido, o “inestable”, así que se debe pivotear • Se computan operaciones BLAS 1 o BLAS 2, pero sabemos que las BLAS 3 (mult. de matrices) es más rápida (clase anterior…) for i = 1 to n-1 A(i+1:n,i) = A(i+1:n,i) / A(i,i) … BLAS 1 (escala un vector) A(i+1:n,i+1:n) = A(i+1:n , i+1:n ) … BLAS 2 (update en rango 1) - A(i+1:n , i) * A(i , i+1:n) Pico BLAS 3 BLAS 2 BLAS 1

  21. Convirtiendo BLAS2 a BLAS3 en GEPP • Bloques • Usado para optimizar mult.matrices • Más difícil en GE por las dependencias de data en GEPP • GRAN IDEA: Updates retrasados • Salvar varios updates consecutivos BLAS2 a una “matriz relegada” • Aplicar muchos updates simultáneamente en una operación BLAS2 • La misma idea funciona en otras partes del alg. lineal • Quedan preguntas por resolver…. • Primer enfoque: necesito escoger un tamaño de bloque b • Algoritmo guarda y aplica b updates • b debe ser suficientemente pequeño de tal forma que la sub matriz activa consistente en b columnas de A quepa en el cache • b debe ser suficientemente grande para que BLAS3 sea rápido

  22. GEPP en bloques (www.netlib.org/lapack/single/sgetrf.f) for ib = 1 to n-1 step b … Procesa b columnas end = ib + b-1 … Apunta al final del bloque de b columns aplica versión BLAS2 de GEPP para obtener A(ib:n , ib:end) = P’ * L’ * U’ … LL denota la parte triangular inferior de A(ib:end , ib:end) + I A(ib:end , end+1:n) = LL-1 * A(ib:end , end+1:n)… update siguientes b filas de U A(end+1:n , end+1:n ) = A(end+1:n , end+1:n ) - A(end+1:n , ib:end) * A(ib:end , end+1:n) … aplica updates retrasados con una mult.matrices … de dimensiones w x b , b x w

  23. Eficiencia del GEPP en bloques

  24. Vista general de LAPACK y ScaLAPACK • Librería estándar para algebra lineal densa y en bandas • Sistemas lineales: A*x=b • Problemas de mínimos cuadrados: minx || A*x-b ||2 • Problemas de Eigenvalues: Ax = lx, Ax = lBx • Descomposición en valores singulares (SVD): A = USVT • Alg. reorganizados para usar lo más posible BLAS3 • Base de librerías matemáticas en muchos sistemas, Matlab, etc … • Aún faltan más mejoras por implementar • Proyectos para maestrías, doctorados?

  25. Performance de LAPACK (n=1000) Performance de Eigen-values, SVD, etc.

  26. Performance de LAPACK (n=100) Eficiencia es mucho menor para matrices más chicas

  27. Paralelizando la GE • Pasos de la paralelización • Descomposición: identificar suficiente trabajo paralelo, pero no demasiado • Asignación: balance de carga entre procesos • Coordinación: comunicación y sincronización • Mapeo: qué procesadores ejecutan qué tareas? • Descomposición • En la parte BLAS 2, casi cada flop en el loop interno puede hacerse en paralelo, así que con n2 procs., necesito 3n pasos paralelos • Esto es grano muy fino, mejor es llamar a mult.matrices locales • Necesario usar mult.matrices paralelo • Asignación • Que procesadores son responsables de qué submatrices? for i = 1 to n-1 A(i+1:n,i) = A(i+1:n,i) / A(i,i) … BLAS 1 (escala un vector) A(i+1:n,i+1:n) = A(i+1:n , i+1:n ) … BLAS 2 (update rango 1) - A(i+1:n , i) * A(i , i+1:n)

  28. Diferentes particiones de data para GE paralelo Mal balance de carga: P0 ocioso luego de n/4 pasos Mejor balance, pero difícil usar BLAS2 o BLAS3 1) Bloque de columnas 1D 2) Columnas 1D cíclicas Se puede negociar balance de carga y performance BLAS2/3 cambiando b, pero la factorización del bloque de columnas es el cuello de botella Difícil de direccionar b 4) Bloques de diagonales 3) Bloques cíclicos de columnas 1D Mal balance de carga: P0 ocioso luego de n/2 pasos El mejor! 6) Bloques cíclicos de filas y columnas 2D 5) Bloques de filas y columnas 2D

  29. Revisión: GEPP con BLAS 3 (con bloques) for ib = 1 to n-1 step b … Procesa matriz de b columnas end = ib + b-1 … Apunta alfinal del bloque de b columnas aplica ver.BLAS2 de GEPP para onbtener A(ib:n , ib:end) = P’ * L’ * U’ … LL denota parte triangular inferior de A(ib:end , ib:end) + I A(ib:end , end+1:n) = LL-1 * A(ib:end , end+1:n)… update sgtes.b filas de U A(end+1:n , end+1:n ) = A(end+1:n , end+1:n ) - A(end+1:n , ib:end) * A(ib:end , end+1:n) … aplicar updates retrasados con una mult.matrices … de dimensiones w x b, x w BLAS 3

  30. Bloque cíclico de Filas y Columnas • Procesadores y bloques de matrices están distribuidos en un array 2d • Array de procesadoresprow x pcol • Bloques de matrices brow x bcol • Paralelismo “pcol” veces en cualquier columna, y llamadas a BLAS2 y BLAS3 en matrices de tamaño brow x bcol • Cuello de botella serial es menos problemático • prow  pcol y brow  bcol posible, es más deseable bcol brow

  31. GE distribuido con Bloques cíclicos 2D • Tamaño de bloque b en el algoritmo y los tamaños de bloque brow y bcol en la partición satisfacen b=bcol. • Regiones sombreadas indican procesadores ocupados con computación o comunicación. • No es necesario tener barrera entre pasos sucesivos del algoritmo, ej: pasos 9, 10, y 11 se pueden hacer en “paralelo”

  32. Distributed GE with a 2D Block Cyclic Layout Contador de columnas Contador de filas

  33. –1 Mult. matrices de verde = verde -azul * rosa

  34. Revisión de MatMul paralela • Se quiere problemas grandes por procesador • PDGEMM = PBLAS mat.mul. • Observaciones: • Para N fijo, si P crece, Mflops crece pero menos que 100% de eficiencia • Si P es fijo, si N crece, Mflops (eficiencia) crece • DGEMM = rutina BLAS • para mat.mul. • Máxima velocidad para PDGEMM • = # Procs * velocidad DGEMM • Observaciones: • Eficiencia siempre es por lo menos 48% • Para N fijo, si P crece, eficiencia disminuye • Para P fijo, si N crece, eficiencia crece

  35. PDGESV = LU paralelo de ScaLAPACK • Como no puede correr más rápido que el loop interno (PDGEMM), hacemos: • Eficiencia = • Veloc(PDGESV)/Veloc(PDGEMM) • Observaciones: • Eficiencia muy por encima de 50% para problemas suficientemente grandes • Para N fijo, si P crece, eficiencia disminuye (igual que en PDGEMM) • Para P fijo, si N crece, eficiencia crece (igual que en PDGEMM) • Según la tabla inferior, el coso de resolver • Ax=b es aprox. la mitad del mat.mul. para matrices suficientemente grandes. • Del contador de flops se deduce que debería ser (2*n3)/(2/3*n3) = 3 veces más rápida, pero la comunicación lo hace algo más lento.

  36. Modelos de performance ScaLAPACK (1) Contador de operaciones de ScaLAPACK tf= 1 tm = a tv = b NB = brow=bcol P = prow = pcol

  37. Modelos de performance ScaLAPACK (2) Comparar Predicciones y Mediciones (LU) (Cholesky)

  38. Escalabilidad de LAPACK y ScaLAPACK • Problemas “unilaterales” son escalables • En GE, A factorizada como el producto de 2 matrices A = LU por medio de pre-multiplicar A por una secuencia de matrices más simples • Asintóticamente es 100% BLAS3 • LU (“Linpack Benchmark”) • Cholesky, QR • Problemas “Bilaterales” som menos escalables • A factorizada como producto de 3 matrices por medio de pre y post multiplicación • 50% BLAS2, no llega a 100% BLAS3 • Eigenproblemas, SVD descomposición de valor singular • Eigenproblemas asimétricos son peores • Problemas de banda estrecha son los peores (para llevarlo a BLAS3 o paralelizar) • Ecuaciones lineales y eigenproblemas • www.netlib.org/{lapack,scalapack}

  39. Algoritmos Recursivos • También usa updates retrasados, pero organizados en forma diferente • Puede explotar particiones de data recursiva • 3x speedups en mínimos cuadrados para matrices altas y delgadas • En teoría, la performance de la jerarquía de memoria es óptima • Referencias • “Recursive Block Algorithms and Hybrid Data Structures,” Elmroth, Gustavson, Jonsson, Kagstrom, SIAM Review, 2004 • http://www.cs.umu.se/research/parallel/recursion/

  40. L A12 A21 A22 GE usando un Algoritmo Recursivo F. Gustavson y S. Toledo Algoritmo LU : 1: Partir la matriz en 2 rectangulos (m x n/2) si solo queda 1 columna, escale por el reciproco del pivot y “return” 2: Aplicar el algoritmo LU a la parte izquierda 3: Aplique transformaciones a la parte derecha (solución triangular A12 = L-1A12 y mult. matrices A22=A22 -A21*A12 ) 4: Aplicar algoritmo LU a la parte derecha La mayor parte del computo es multiplicar matrices de tamaño n/2, n/4, n/8, … Source: Jack Dongarra

  41. Factorizaciones recursivas • Tan precisas como los métodos convencionales • Mismo número de operaciones • Se forman bloques de tamaño variable en forma automática • Sólo operaciones BLAS nivel 1 y 3! • Simplicidad de expresión • Potencialmente eficientes pero a la vez “indiferente” al caché • Pero no se debería llevar la recursión hasta el extremo de 1 columna! • La formulación recursiva es sólo un reajuste de algoritmos más viejos del LINPACK

  42. Recursive LU Dual-processor LAPACK Recursive LU LAPACK Uniprocessor Source: Jack Dongarra

  43. Limites de los algortimos recursivos • Dos clases de descomposiciones de matrices densas • “Unilaterales” • Secuencia de operaciones sencillas aplicadas en la izquierda de la matriz • GE: A = L*U o A = P*L*U • GE simétrica: A = L*D*LT • Cholesky, sólo matrices simétricas: A = L*LT • Descomposición QR para mínimos cuadrados: A = Q*R • Puede acercarse a 100% BLAS 3 • Susceptible de recursión • “Bilaterales” • Secuencia de operaciones sencillas aplicadas en ambos lados, alternadamente • Algoritmos de Eigenvalores, SVD • Por lo menos ~25% BLAS 2 • Enfoque recursivo imposible, difícil? • Algunos progresos la última decada: SVD (25% vs 50% BLAS2)

  44. Algoritmos “Out of Core” Out-of-core significa que la matriz reside en disco; muy grande para entrar completa en la RAM Mucho más difícil evitar efecto de latencia del disco QR más fácil que LU por que no hace falta pivotear para QR Source: Jack Dongarra

  45. Resolviendo Sist.Ec.Linealesprovenientes de Ec.Dif.Parc.

  46. Repaso de “Fuentes de paralelismo” • Sistemas de eventos discretos: • Ejs: “Game of Life,” simulación de circuitos lógicos. • Sistemas de partículas: • Ejs: bolas de billar, simulación de dispositivos semiconductores, galaxias. • Variables agregadas dependientes de parámetros continuos: • EDOs, ejs: simulación de circuitos (Spice), mecánica estructral, cinética química. • Variables contínuas dependientes de parámetros continuos: • EDPs, ejs: calor, elasticidad, electrostática. • Un fenómeno dado puede ser modelado en múltiples niveles. • Muchas simulaciones combinan más de una técnica.

  47. Repaso: Resolviendo EDPs ° Problemas hiperbólicos (ondas): • Onda sonora (posición, time) • Usan iteración en el tiempo en forma explícita • Solución en cada punto depende de los vecinos en Δt anterior • Problemas elípticos (estado estable): • Potencial electrostático (posición) • Todo valor depende de todos lo demás valores • Localidad es más difícil de hallar que en los problemas hiperbólicos • Problemas parabólicos (dependientes del tiempo): • Temperatura (posición, tiempo) • Involucra una solución elíptica cada Δt • Nos enfocaremos en los elípticos • Ejemplo canónico es la ec. Poisson 2u/x2 + 2u/y2 + 2u/z2 = f(x,y,z)

  48. Ec. Poisson aparece en muchos lados 3D: 2u/x2 + 2u/y2 + 2u/z2 = f(x,y,z) f representa las fuentes; además se necesita condiciones de frontera 2D: 2u/x2 + 2u/y2 = f(x,y) 1D: d2u/dx2 = f(x) • Potencial electrostático o gravitacional: Potencial (posición) • Flujo de calor: Temperatura (posición, tiempo) • Difusión: Concentración (posición, tiempo) • Flujo de fluidos: Velocidad, Presión ,Densidad (posición, tiempo) • Elasticidad: Tensión, estiramiento (posición, tiempo) • Variaciones de Poisson tienen coeficientes variables

  49. Relación de ec. Poisson con la Gravedad, Electrostática • Fuerza en la partícula en (x,y,z) debido a una partícula en (0,0,0) • 3D: -(x,y,z)/r3, where r = (x2 +y2 +z2 ) • 2D: -(x,y)/r2, where r = (x2 + y2 ) • Fuerza es gradiente del potencial V (con signo opuesto) • 3D: V = -1/r, V = (V/x, V/y, V/z) • 2D: V = log r, V = (V/x, V/y) • V satisface la ec. de Poisson (tarea)

  50. Ec. Poisson en 1D: 2u/x2 = f(x) Discretizar d2u/dx2 = f(x) en la malla regular ui = u(i*h) Para conseguir [ u i+1 – 2*u i + u i-1 ] / h2 = f(x) Se resuelve como Tu = -h2 * f Para array desconocido u donde 2 -1 -1 2 -1 -1 2 -1 -1 2 -1 -1 2 Grafo y “estencil” -1 2 -1 T =

More Related