Ordenação Lexicográfica: Aumento de Complexidade versus Redução de Recurso Computacional

 Durante a implementação de um código para modelar a transferência de calor por radiação em um domínio bidimensional, me deparei com o seguinte problema: mesmo levando em conta as simplificações, a grandeza fundamental do problema, a intensidade da radiação, precisaria ser descrita por quatro dimensões, ou seja, é uma variável quadridimensional. As dimensões são: duas coordenadas espaciais: r e z (o problema é em coordenadas cilíndricas) e duas coordenadas que definem a direção de propagação da radiação: o ângulo polar theta e o ângulo azimutal phi.

Não é incomum que o modelo numérico de problemas de CFD e CHT produzam matrizes esparsas, ou seja, matrizes onde apenas a diagonal principal e as diagonais adjacentes são formadas por elementos não nulos. Neste caso, o uso de matrizes para computar a relação entre um ente discreto (ponto, volume ou elemento) e os vizinhos faz com que seja necessário alocar muita memória que não será usada. Além de limitar o tamanho do problema que o hardware resolver, isso aumenta o tempo de processamento sobremaneira. Quando eu estudei mecânica dos fluidos computacional, aprendi uma técnica simples e vantajosa para reduzir a memória requerida na solução desses problemas cujo modelo numérico recai em matrizes esparsas: armazenar os coeficientes que relacionam os entes discretos em vetores e programar o código de forma a buscar os coeficientes nos vetores que correspondem à respectiva posição na 'matriz original'. Para isso ser feito em problemas multidimensionais, onde a matriz dos coeficientes possui várias 'diagonais' não nulas (diagonais paralelas à diagonal principal da matriz), é necessário empregar a ordenação lexicográfica.

A ordenação lexicográfica é uma indexação da posição dos elementos de uma matriz (ou ente de ordem superior) de uma forma previamente definida em um vetor (ou ente de ordem inferior ao original). Toda vez que for necessário varrer a matriz em busca de determinada posição será necessário correr a ordenação lexicográfica predeterminada. É um pequeno prejuízo em termos de processamento, mas que pode resultar em grande economia de memória.

Exemplo: Temos uma matriz no centro da Figura 1. Nela, cada posição é denotada pelo par ordenado (i,j) que identifica a posição na matriz conforme a orientação dos eixos coordenados (x,y). O número de elementos na direção x é Ni e na direção y é Nj. Ao redor da matriz estão mostradas algumas opções de ordenação lexicográfica com as respectivas fórmulas de ordenação, onde r é o respectivo índice da posição da matriz em termos do vetor já ordenado. Note que a orientação da ordenação é imaterial, ou seja, independe da direção.


Figura 1 - Exemplos de ordenação lexicográfica de malhas bidimensionais uniformes


Agora chegamos ao ponto interessante desta postagem! Alguns meses atrás eu iniciei a programação da versão bidimensional do Método das Ordenadas Discretas (DOM do termo em língua inglesa Discrete Ordinates Method). Este método é usado para calcular a transferência de calor por radiação em meios participantes. Minha experiência até então era na programação da versão unidimensional, mas as coisas complicaram ao me deparar com a necessidade de implementar a versão bidimensional em coordenadas cilíndricas.

Embora não seja uma postagem para explicar o procedimento do DOM, é importante ressaltar que ele requer uma discretização do espaço (similar aos demais métodos empregados em CFD) e uma discretização angular (em cada elemento de volume é necessário resolver a equação da transferência radiativa em várias direções).

Na versão bidimensional as direções ordenadas são mostradas conforme a Figura 2 a seguir, onde csi e mi são os cossenos diretores das direções ordenadas representadas unindo-se a origem com cada um dos 24 pontos mostrados, como feito para a direção indicada pela seta vermelha. 

Figura 2 - Representação da Aproximação angular S6 para problemas bidimensionais

O número de direções ND é calculado por ND=2^D*SN*(SN+2)/8, onde D={1,2,3} é a dimensão do problema e SN é a ordem da aproximação angular, que neste exemplo é 6, pois há seis níveis do cosseno diretor principal (csi) constante, conforme mostrado. Desta forma, a Figura 2 mostra a representação da aproximação S6. Em versão tridimensional, esta e outras aproximações são mostradas na Figura 3 a seguir.

Figura 3 - Octante de uma esfera de raio unitário mostrando as direções positivas de quadraturas tridimensionais simétricas e rotacionalmente invariantes.

No DOM, a varredura do domínio ocorre em dois laços de programação do tipo "do / end do": o mais externo varre as direções e o mais interno as posições no espaço. Desta forma, dependendo da direção, o laço interno visita os elementos de volume seguindo orientações lexicográficas como aquelas exemplificadas na Figura 1. Porém o próprio laço externo, que corre as direções, não possui um número de direções constante em cada nível, conforme se pode ver na Figura 2. Nota-se que para csi = 1 há duas direções, enquanto para csi = 2 há 4, etc. Entretanto, esse fato não impede que seja encontrada uma lógica que permita correr as direções mesmo que seu número varie de acordo com o nível de csi. Para explicar essa lógica, podem ser empregados dois contadores auxiliares i e j, dispostos conforme a Figura 4 a seguir. 

Figura 4 - Proposição dos índices auxiliares i e j, usados para correr a ordenação lexicográfica no DOM.

A ordenação lexicográfica desejada é mostrada na Figura 5, onde p = {1, 2, ... ND} é o índice que corre os níveis de csi e q = {1, 2, ..., Np} é o índice que corre os níveis na direção mi, onde Np é o número de direções do nível p. Ela inicia no primeiro nível de csi e menor índice possível de mi e segue até atingir o maior índice de mi. Ao atingi-lo, o método incrementa o índice de csi e o processo continua conforme mostrado pelas setas azuis. Na Figura 5 as siglas Q1, ..., Q4 representam os quadrantes 1 até 4, usados apenas para ajudar a situar onde devem ser inseridos os cálculos do método DOM entre os trechos de código da ordenação lexicográfica .

Figura 5 - Ordenação lexicográfica para correr as direções ordenadas.

Usando a linguagem FORTRAN95, uma das formas de programar a ordenação lexicográfica desejada é usando o seguinte trecho de código cada vez que se faz necessário correr as direções ordenadas:

i = -SN/2 - 1

pto1 = 0

r = 0

do p = 1, SN/2

    pto = pto1

    i = i + 1

    j = SN/2 + 1 - abs(i)

    j = sign(j,i)

    do q = 1, p

        pto = pto + 1

        r = r + 1

       

        ! INSERIR CÁLCULOS DO QUADRANTE 3 AQUI

       

        j = j + 1

    end do

    pto1 = pto

    j = j + 1

    do q = p+1, 2*p

        r = r + 1

       

        INSERIR CÁLCULOS DO QUADRANTE 4 AQUI

       

        pto = pto - 1

        j = j + 1

    end do

end do

pto1 = SN(SN+2)/8

pto2 = 0

i = i + 1

do p = SN/2+1, SN

    pto2 = pto2 + 1

    pto = pto1 - (SN/2 - pto2 + 1)

    i = i + 1

    j = SN/2 + 1 - abs(i)

    j = - sign(j,i)

    do q = 1, SN-p+1

        pto = pto + 1

        r = r + 1

       

        ! INSERIR CÁLCULOS DO QUADRANTE 2 AQUI

       

        j = j + 1

    end do

    j = j + 1

    do q = SN-p+2, 2*(SN-p+1)

        r = r + 1

       

        INSERIR CÁLCULOS DO QUADRANTE 1 AQUI

       

        pto = pto - 1

        j = j + 1

    end do

    pto1 = pto

end do


onde o sinal = representa atribuição (e não igualdade) e os índices p, qr, i, j, pto e pto1 são contadores (variáveis tipo inteiro) e SN é a ordem da aproximação angular (também do tipo inteiro). A função abs(i) toma apenas o valor da variável i em módulo e a função sign(i,j) aplica na variável i o sinal do número contido na variável j.

A armazenagem das variáveis multidimensionais como a intensidade é feita usando o índice r, que é uma combinação dos índices p e q. No exemplo acima, o problema abordado é bidimensional, portanto = {1, 2, ..., 24}. Como a intensidade é uma variável quadridimensional, ao ser aplicada a ordenação lexicográfica na discretização angular, esta passa a ser tridimensional e se também for aplicada na discretização espacial, então a intensidade passa a ser uma variável bidimensional apenas. Ou seja, são reduzidas duas dimensões no modelo numérico!

É interessante notar que apesar de o código ficar um pouco mais difícil do programador interpretar, o uso da ordenação lexicográfica permite elaborar lógicas relativamente complexas de forma relativamente fácil. Um exemplo é a aplicação da condição de contorno de reflexão especular, onde a intensidade atingindo uma superfície em uma dada direção deve ser refletiva em outra direção, que forma o mesmo ângulo em relação à normal à superfície. Uma matriz ND x 2 ou ND x 3 já dá conta de armazenar a informação de qual direção reflete em qual outra direção, ocupando pouca memória.

Entretanto, a maior justificativa para o uso da ordenação lexicográfica é o fato que requer menos memória para armazenar os cossenos diretores (3 componentes) e pesos (1 valor) das direções. No exemplo supracitado, onde SN = 6, a ordenação sem ordem lexicográfica usará uma matriz com 6 x 6 = 36 linhas e 4 colunas (onde várias posições não serão usadas), enquanto que com a ordenação lexicográfica aqui proposta a mesma matriz usará 24 linhas por 4 colunas, que é a quantidade mínima possível. É uma redução de 33,33% e quanto maior a ordem da aproximação, mais significativa a quantidade de memória economizada com essa técnica. No caso SN = 100, a redução é de 51%. A Figura 6 mostra a razão entre a memória ocupada pela ordenação lexicográfica e a ordenação convencional para um problema bidimensional onde o DOM é empregado.

Figura 6 - Razão entre a memória ocupada pela ordenação lexicográfica e a ordenação convencional em um problema bidimensional de radiação térmica usando o DOM

Comentários

Postagens mais visitadas deste blog

Engenharia Reversa do "Boundary Layer Module" - Parte 1

Engenharia Reversa do "Boundary Layer Module" - Parte 2