Neste trabalho aplica‐se o método sem malha natural neighbour radial point interpolation method (NNRPIM) à análise unidimensional de vigas laminadas, considerando a teoria de Timoshenko.
O NNRPIM combina o conceito matemático dos vizinhos naturais com a interpolação radial pontual. Os diagramas de Voronoï permitem impor a conectividade nodal e construir a malha de fundo para efeitos de integração, por intermédio das células de influência. É apresentada a construção das funções de interpolação NNRPIM, sendo, para estas, usada a função de base radial multiquadrática. As funções de interpolação geradas possuem continuidade infinita e a propriedade de delta Kronecker, o que facilita a imposição das condições de fronteira, uma vez que estas podem ser impostas com o método da imposição direta, tal como no método dos elementos finitos (FEM).
De modo a obter o campo de deslocamentos e de deformações, a teoria de deformação de Timoshenko para vigas sujeitas a esforços transversos é considerada. Vários exemplos numéricos de vigas isotrópicas e vigas laminadas são apresentados de modo a demonstrar a convergência e a exatidão da aplicação proposta. Os resultados obtidos são comparados com soluções analíticas disponíveis na literatura.
In this work, a meshless method, “natural neighbour radial point interpolation method” (NNRPIM), is applied to the one‐dimensional analysis of laminated beams, considering the theory of Timoshenko.
The NNRPIM combines the mathematical concept of natural neighbours with the radial point interpolation. Voronoï diagrams allows to impose the nodal connectivity and the construction of a background mesh for integration purposes, via influence cells. The construction of the NNRPIM interpolation functions is shown, and, for this, it is used the multiquadratic radial basis function. The generated interpolation functions possess infinite continuity and the delta Kronecker property, which facilitates the enforcement of boundary conditions, since these can be directly imposed, as in the finite element method (FEM).
In order to obtain the displacements and the deformation fields, it is considered the Timoshenko theory for beams under transverse efforts. Several numerical examples of isotropic beams and laminated beams are presented in order to demonstrate the convergence and accuracy of the proposed application. The results obtained are compared with analytical solutions available in the literature.
Recentemente, os métodos sem malha [1] têm‐se afirmado como alternativa numérica ao método dos elementos finitos (MEF)[2], permitindo a resolução de problemas com geometrias muito complexas, apresentando um menor tempo de processamento [1,3]. Distinto do MEF, onde o domínio é discretizado numa malha de elementos, nos métodos sem malha o domínio do problema é discretizado numa malha de nós e a conectividade nodal é obtida por domínios de influência. Estes domínios, ao contrário dos elementos no MEF, podem e devem sobrepor‐se. É esta sobreposição e a escolha de domínios de influência convenientes que permitem aos métodos sem malha, mesmo perante malhas nodais muito irregulares, obter bons resultados.
Numa primeira fase foram desenvolvidos vários métodos sem malha com recurso a funções aproximadoras –como o element free galerkin method (EFGM) [1] ou o meshless local Petrov‐Galerkin method[4]– no entanto, apesar destes terem sido aplicados com sucesso à mecânica computacional, alguns problemas ficaram por resolver, tal como a imposição das condições de fronteira, devido à falta da propriedade de delta Kronecker. Esta é uma das consequências diretas do uso de funções aproximadoras em vez de funções interpoladoras. Assim, de forma a dar resposta a esta limitação, novos métodos sem malha foram desenvolvidos na última década, entre os quais o point interpolation method (PIM) [5,6], o natural neighbour finite element method (NNFEM) [7,8] e o meshless finite element method (MFEM) [9], fazendo uso de funções interpoladoras.
O PIM é um método bastante atrativo pelo uso de funções interpoladoras, gerando assim funções de forma com a propriedade de delta Kronecker. A construção destas funções de forma, tal como as suas derivadas parciais, é um processo simples, quando comparado com o EFGM. Entretanto o PIM evoluiu e em vez de utilizar funções de base polinomiais foi proposto, e aplicado com sucesso, o uso de funções de base radial [10,11], sendo criado o radial point interpolation method (RPIM). Nos primeiros trabalhos sobre o RPIM as funções de base radial utilizadas foram a multiquadrática e a gaussiana [10,11]. Este utiliza o conceito de «domínio local de influência» em vez de «domínio global de influência», gerando assim matrizes de rigidez esparsas e em banda, mais ajustadas a problemas de geometrias complexas.
O natural element method (NEM) é atualmente um dos métodos sem malha mais populares, isto porque as funções de forma são construídas com base nos vizinhos naturais, através de funções interpoladoras. A construção destas funções interpoladoras baseia‐se em conceitos matemáticos, como são o caso dos diagramas de Voronoï [12].
A necessidade de uma malha de fundo para efeitos de integração é também um tema de grande discussão nos métodos sem malha. O método apresentado neste trabalho dá uma resposta original e eficiente a essa necessidade.
Neste trabalho o natural neighbour radial point interpolation method (NNRPIM) [13] é aplicado à resolução de problemas de vigas, considerando a teoria da deformação de Timoshenko (teoria de deformação de corte de primeira ordem). Recorrendo às células de Voronoï, um conjunto de células de influência é criado a partir de uma discretização nodal. Para efeitos de integração é criada uma malha de fundo com pontos de quadratura distribuídos segundo a regra gaussiana.
O NNRPIM utiliza funções de base radial (RBF) para construir as funções de interpolação. Na literatura é possível encontrar vários trabalhos sobre a obtenção dos parâmetros de forma ótimos das RBF [14–16] que conduzem à melhor solução. Relativamente ao NNRPIM, é possível encontrar na literatura estudos sobre os parâmetros de forma ótimos para o problema da mecânica dos sólidos bidimensional e tridimensional [13], assim como para a formulação de Reissner‐Mindlin para placas e laminados [17]. Os parâmetros de forma obtidos no trabalho de Dinis et al. [13,17] diferem dos valores obtidos no trabalho de Wang and Liu [11], no entanto, os parâmetros de forma obtidos por Dinis et al. [13,17] conduzem a soluções mais próximas da solução considerada como exata.
É notório, após consultar a literatura disponível, o interesse contínuo por parte da comunidade académica no sentido de desenvolver modelos matemáticos eficientes na previsão da resposta estática e dinâmica deste tipo de estruturas unidimensionais. É possível consultar em [18] uma revisão detalhada sobre algumas das teorias de deformação de vigas e placas usadas na análise estática e dinâmica. O problema de vigas espessas encontra‐se bem estudado no campo dos elementos finitos [19,20]. No entanto, no domínio dos métodos sem malha este é ainda um tema pouco explorado [21–24].
2Vizinhos naturaisSurgido em 1980, por Sibson [25], o conceito de vizinhos naturais permite impor uma conectividade nodal celular [26], traduzindo‐se no conjunto de nós mais próximos do nó em análise que, no caso unidimensional representado na figura 1(a), serão o nó mais à esquerda e o nó mais à direita.
Os vizinhos naturais são determinados recorrendo ao diagrama de Voronoï, sendo que este é constituído por um conjunto de células de Voronoï. A cada uma dessas células está associado apenas um nó, para além do conjunto de todos os pontos que estão mais próximos desse nó do que de qualquer outro.
Na figura 1 são exemplificadas as diferentes etapas para a determinação do diagrama de Voronoï. Assim, numa primeira fase, o domínio do problema é discretizado numa malha nodal com N nós, N={n1, n2, n3, ..., nN}.
Posteriormente, é construída uma célula provisória, VI*, que contem os nós mais próximos do nó em estudo (nI), ou seja, que contém os seus nós vizinhos (nI‐1 e nI+1), figura 1(b).
De seguida, é determinada a célula de Voronoï, figura 1 (c). Esta tem por base a célula provisória VI*, no entanto, a distância entre os seus limites e nI será obtida por ε(xJ,xI)2=dIJ=(xI−xJ)22, sendo ¿ a norma euclidiana de xJ e xI, com J≠I, onde J é um nó vizinho de I.
No final, obtém‐se o diagrama de Voronoï do domínio nodal [13,27], figura 1(d).
O diagrama de Voronoï permite impor a conectividade nodal e construir a malha de fundo para efeitos de integração [27].
2.1Células de influência e conectividade nodalNo caso do NNRPIM as células de influência servem para impor a conectividade nodal e permitem a criação da malha de integração, sendo constituídas por um conjunto de nós na vizinhança do ponto de interesse [13,26,27]. Três tipos distintos de células de influência são seguidamente apresentados na figura 2 para o nó I (nI).
Como é visível, as células de influência de terceiro grau, figura 2(c), apresentam uma maior conectividade nodal. Espera‐se, assim, com as células de influência de terceiro grau, obter resultados mais próximos da solução analítica.
2.2Integração Numérica no natural neighbour radial point interpolation methodO esquema de integração usado no NNRPIM 1D consiste, numa primeira fase, em definir as células de Voronoï, figura 3(a). Posteriormente, são distribuídos na partição celular definida pela distância dIJ o número pretendido de pontos de quadratura, figura 3(b), que permitirão a integração numérica das funções NNRPIM. Essa distribuição é feita segundo a regra de Gauss, onde cada ponto de integração tem um peso associado.
No caso da integração nodal, os pontos de integração são os próprios nós, pelo que essa apresentará um aspeto semelhante à figura 3(a). Nesse caso, o peso correspondente ao ponto de integração é a distância l.
Na secção 5 é apresentado um estudo onde é obtido o número de pontos de integração por partição celular, necessários e suficientes para integrar numericamente o problema da viga de Timoshenko.
3Funções de forma natural neighbour radial point interpolation methodNeste trabalho foram usadas funções interpolantes, do tipo pontual radial – RPI [10]. Considere‐se uma função u(x), definida num domínio Ω e discretizada por um conjunto de N nós, onde se assume que apenas os nós pertencentes à célula de influência afetam a função u(x), e que esta passa por todos os nós, usando uma função de base radial. Assim, e considerando o valor da função num ponto de interesse xI, tem‐se que
Sendo que RixI é a função de base radial (RBF), n é o número de nós no interior da célula de influência de xI, e ai(xI) e bj(xI) são coeficientes não constantes de RixI e pjxI, respetivamente. O monómio j da base polinomial é definido por pjxI e m é o número de monómios da base. Os vetores da equação (1) definem‐se por,
Nas funções de base radial (RBF) a variável é a norma euclidiana (rij) entre o ponto de interesse xI e o nó vizinho, rij=(xI−xi)2.
Neste trabalho é usada a função de base radial multiquadrática, R(rij)=(rij2+c2)p.
Esta foi proposta por Hardy [27,28], e os parâmetros c e p são os parâmetros de forma que carecem de otimização, de modo a maximizar a performance do método. Na secção 5 é apresentado um estudo da otimização dos parâmetros c e p.
Por forma a garantir que a matriz de interpolação das funções de base radial é invertível, as bases polinomiais usadas têm de ser completas [13]. Assim, neste trabalho serão usadas as seguintes bases polinomiais 1D:
Para que a solução obtida seja única, a base polinomial tem ainda de satisfazer os requerimentos extra, [13,27]
Assim, uma nova equação matricial pode ser escrita
e desenvolvida,
As funções teste são obtidas substituindo a equação (9) na equação (1),
Sendo φ(x) a função de interpolação definida por
A figura 4 mostra esquematicamente as funções de interpolação para uma célula de influência de primeiro, segundo e terceiro grau, respetivamente.
De notar que as funções de interpolação apresentam valor máximo para o nó em estudo e que passam obrigatoriamente por todos os restantes nós pertencentes à célula de influência.
4Equações de equilíbrioNeste trabalho, a teoria de deformação de primeira ordem é aplicada a vigas (fig. 5). A teoria de Timoshenko, teoria em que este trabalho se baseia, assenta em 3 hipóteses, que são:
- •
Os deslocamentos segundo z de todos os pontos da secção transversal da viga são pequenos e iguais aos deslocamentos do eixo da viga;
- •
A deformação segundo y é nula;
- •
As secções transversais retas normais ao eixo antes da deformação da viga permanecem planas após a deformação, mas não necessariamente normais ao eixo da viga.
Aplicação das hipóteses da teoria de Timoshenko a uma viga laminada30.
Os resultados obtidos por estas hipóteses são, no entanto, uma melhor aproximação quando se tratam de vigas com uma maior relação espessura/vão, uma vez que esta teoria de deformação permite considerar a deformação de corte.
O campo de deslocamentos pode ser obtido a partir da figura 5. Assim,
No caso da teoria de Timoshenko, o vetor das deformações reduz‐se a:
Assim,
O que simplificando,
As deformações podem ainda ser representadas como,
Onde L é um operador diferencial e u é o vetor campo de deslocamentos. Assim,
O valor da tensão na camada i é obtido por aplicação da lei de Hooke,
Onde ci é a matriz constitutiva do material da camada i, definida por,
Sendo Ei o módulo de Young da camada i e νi o valor do coeficiente de Poisson respetivo, k é o fator de correção ao corte de valor igual a 5/6. Considere o sólido de domínio Ω e fronteira Г. As equações de equilíbrio são expressas por
Sendo ∇ o gradiente, σ a tensão e b o conjunto das forças externas aplicadas no corpo. Aplicando a forma fraca de Galerkin,
Sendo q a força de corpo, e t a força externa aplicada na fronteira natural Γ.
Substituindo a equação (16) no primeiro termo da equação (21) obtém‐se,
Considerando a equação (10),
Sendo B a matriz de deformabilidade, a qual pode ser decomposta,
Onde Bf, BA e BC são as componentes de B correspondentes à flexão, ao esforço axial e ao corte, respetivamente.
No problema em estudo as matrizes de deformabilidade correspondentes ao esforço axial e ao corte podem ser agrupadas,
Assim, a matriz de rigidez K do problema da viga de Timoshenko pode ser obtido a partir de,
Sabendo que dΩ=dA·dz e dA=b·dx, uma vez que a largura da viga é constante, a equação (26) pode ser reescrita da seguinte forma,
Considerando a propriedade distributiva e a largura da viga unitária (b=1),
Sendo k o número de camadas do laminado e zi e zi+1 respetivamente as coordenadas das faces inferior e superior da camada i, figura 5.
O segundo termo da equação (21) pode ser desenvolvido considerando a equação (10),
sendo a matriz H definida por,
e q por,onde N é o esforço axial, V o esforço de corte e M o momento fletor.
Raciocínio idêntico pode ser seguido para o terceiro termo da equação (21),
sendo Fq=∫HqdΩ e Ft=∫HtdΓ. A equação (21) pode então ser apresentada na forma matricial,sendo K a matriz de rigidez global, F=Fq+Ft o vetor global dos momentos e das forças externas aplicadas e U o campo de deslocamentos. Assim, o sistema linear de equações pode ser apresentado na forma matricial, KU=F.
As condições de fronteira são impostas diretamente na matriz de rigidez K, assim como no MEF, uma vez que as funções de interpolação radiais pontuais possuem a propriedade de delta Kronecker [17].
5Otimização dos parâmetros de forma da Radial Basis FunctionNesta secção são determinados os parâmetros de forma p e c ótimos, e o número de pontos de quadratura suficientes para integrar o problema de vigas, para cada um dos tipos de células de influência apresentados na secção 2.1.
Foi analisada uma viga isotrópica (E=1.000Pa, ν=0,3) de dimensões L=10m e h=1m (L/h=10) simplesmente apoiada em ambas as extremidades e sujeita a uma carga uniformemente distribuída q=q0, sendo q0=1N/m. Nesta análise a viga foi discretizada em 33 nós, distribuídos numa malha regular e numa outra irregular.
O erro médio do deslocamento transversal é obtido através de,
onde δw NNRPIMi e δw Analíticai são os deslocamentos transversais obtidos no nó i pelo NNRPIM e pela equação analítica para vigas finas, respetivamente, e N é o número total de nós que discretiza o domínio do problema.
Fez‐se variar os valores do parâmetro c entre 0,0001 e 10,0001 e estudou‐se a variação do erro médio do deslocamento vertical (εw) para a viga em questão (malha regular e malha irregular), considerando diferentes bases polinomiais – ausência de base polinomial (p0), base polinomial constante (p1), linear (p2) e quadrática (p3).
Na figura 6 apresentam‐se os resultados obtidos para células de influência de segundo grau, v2. Resultados semelhantes são obtidos quando células de influência de primeiro grau, v1, e de terceiro grau, v3, são consideradas.
Na tabela 1 encontra‐se um quadro resumo dos parâmetros c ótimos obtidos para cada uma das formulações.
Considerando os parâmetros c da tabela 1, procedeu‐se a uma nova análise do problema para a mesma viga (L/h=10), fazendo‐se variar o parâmetro p entre 0,0001 e 5,0001. Na figura 7 apresentam‐se os resultados, novamente, para células de influência de segundo grau.
Na tabela 2 são apresentados os parâmetros p ótimos obtidos para cada uma das formulações.
Considerando, para cada formulação, os parâmetros c e p da tabelas 1 e 2, respetivamente, procedeu‐se a uma nova análise do problema para a mesma viga, fazendo‐se variar o número de pontos de quadratura entre 0‐10 em cada partição celular. Na figura 8 são apresentados os resultados para células de influência de segundo grau.
Na tabela 3 são apresentados os valores do número de pontos de quadratura suficientes para as diferentes formulações.
Pela análise das tabelas anteriores (tabelas 1–3) é ainda possível observar que não se obtiveram resultados para a formulação v1p3, pois o uso de uma base polinomial quadrática, para o caso do problema unidimensional, exige no mínimo 3 nós no interior de cada célula de influência. Para valores inferiores o problema torna‐se indeterminado, sendo a matriz geométrica G não‐invertível.
Quanto ao valor do parâmetro ótimo c, ficou claro que este deverá tomar o valor 0,0001, independentemente da formulação usada. O mesmo se sucedeu para o número ótimo de pontos de integração, sendo que neste caso o seu valor ótimo será 2.
O valor do parâmetro p ótimo, no entanto, não parece ser tão geral e sofre algumas variações nas diferentes formulações.
Após a análise dos resultados obtidos, verifica‐se que a formulação que apresenta melhores resultados é a formulação v2p0, uma vez que é a que aparenta ser mais estável e ter um menor erro. Assim, daqui em diante esta será a formulação usada para os diferentes problemas de vigas que serão considerados.
6Exemplos numéricos6.1Estudo de vigasNesta secção procede‐se ao estudo da convergência de 4 vigas simplesmente apoiadas, com diferentes tipos de cargas aplicadas, utilizando o NNRPIM. São obtidos os resultados do deslocamento transversal no centro da viga para malhas regulares e irregulares, em função do número de divisões da viga (número de nós=número de divisões+1). As soluções NNRPIM são comparadas com as soluções analíticas. São também obtidos, para cada tipo de viga, com a respetiva carga, o erro médio do deslocamento transversal através da equação (34). Os resultados do estudo de convergência são apresentados nas figuras 10–13.
O estudo foi efetuado para vigas com as dimensões L=10m e h=1m, sendo as propriedades materiais E=1.000Pa e ν=0,3.
Os 4 tipos de carregamentos usados estão representados na figura 9, sendo q0=1N/m.
Os resultados apresentados nas Figuras 10(a), 11(a), 12(a), 13(a) e tabela 4 encontram‐se sob a forma de deslocamentos transversais normalizados. Sendo essa normalização efetuada por,
Valores dos deslocamentos normalizados para diferentes relações L/h
| Deslocamentos normalizados | ||||||
| L/h=5 | L/h=10 | L/h=20 | L/h=50 | L/h=100 | ||
| v2p0 regular | 2,812 | 2,578 | 2,519 | 2,503 | 2,499 | |
| v2p0 irregular | 2,812 | 2,578 | 2,520 | 2,502 | 2,496 | |
| Analítica | 2,500 | 2,500 | 2,500 | 2,500 | 2,500 | |
| v2p0 regular | 17,185 | 16,015 | 15,722 | 15,638 | 15,618 | |
| v2p0 irregular | 17,188 | 16,017 | 15,724 | 15,635 | 15,601 | |
| Analítica | 15,625 | 15,625 | 15,625 | 15,625 | 15,625 | |
| v2p0 regular | 8,593 | 8,008 | 7,861 | 7,819 | 7,809 | |
| v2p0 irregular | 8,594 | 8,008 | 7,862 | 7,818 | 7,800 | |
| Analítica | 7,813 | 7,813 | 7,813 | 7,813 | 7,813 | |
| v2p0 regular | 14,008 | 13,033 | 12,789 | 12,719 | 12,702 | |
| v2p0 irregular | 14,011 | 13,035 | 12,790 | 12,717 | 12,708 | |
| Analítica | 12,708 | 12,708 | 12,708 | 12,708 | 12,708 | |
Na tabela 4 são apresentados os valores dos deslocamentos transversais no centro da viga obtidos considerando diferentes relações L/h e uma malha única de 1.025 nós.
Apresentam‐se de seguida os diagramas de esforços transversos e momentos fletores para as vigas estudadas com dimensões L=10m e h=1m (fig. 9). Foi usada uma malha regular com 1025 nós.
Da análise das figuras 14–17 constata‐se que a distribuição do campo de momentos e esforços transversos é suave ao longo de toda a viga.
De modo a estudar o comportamento do NNRPIM em situações de ocorrência do fenómeno de retenção ao corte, fez‐se variar a altura da viga h entre h=0,0001 e h=5 m, mantendo‐se o comprimento constante, L=10.
Assim, na figura 18 são apresentados os resultados obtidos considerando as condições de fronteira essenciais e naturais indicadas.
Analisou‐se o problema recorrendo a 2 esquemas de integração distintos, a integração nodal, p.q.=0, e a integração utilizando pontos de quadratura de Gauss (secção 2.2), p.q.=k, sendo k o número de pontos de quadratura por partição de célula de integração.
6.2Estudo de laminadosNesta secção procede‐se ao estudo do problema de vigas laminadas. Analisam‐se vigas com condições de fronteira distintas, sujeitas a 2 tipos de carregamentos, considerando diferentes tipos de laminados (fig. 19).
Na figura 20 são apresentadas genericamente as propriedades geométricas e materiais das vigas laminadas em estudo. Sendo Ei o módulo de elasticidade na direção material i e Gij o módulo de distorção entre as direções materiais i e j.
O deslocamento transversal no centro da viga laminada obtido com o NNRPIM, wmax, é normalizado com base na expressão,
onde q0=1N/m. Para o caso da carga uniformemente distribuída q=q0, já para a carga pontual Q=q0·L. Os resultados normalizados obtidos pelo NNRPIM são comparados com a solução analítica [29]. Os resultados obtidos encontram‐se nas figuras 21–23 para as vigas laminadas a 0°, 90° e 0°/90°/90°/0°, respetivamente. Na tabela 5 é apresentada uma comparação entre algumas das soluções obtidas com o NNRPIM e soluções disponíveis na literatura obtidas por outros autores.
Na figura 24 é possível observar os diagramas para a tensão normal e de corte da viga laminada a 0°/90°/90°/0°.
Fez‐se agora variar a espessura da viga h entre h=0,10 e h=5 mantendo‐se o comprimento L=10.
Na figura 25 são apresentados os deslocamentos transversais normalizados obtidos para as diferentes vigas laminadas consideradas, em função das suas condições de fronteira e tipo de carregamento.
7Conclusões
Neste trabalho foi estendido à análise de vigas laminadas um método sem malha, o NNRPIM. Por forma a obter a conectividade nodal é construído o diagrama de Voronoï da malha nodal que discretiza o domínio do problema. Com base nas células de Voronoï obtidas é construída a malha de integração.
Distintas formulações NNRPIM foram estudadas, considerando células de influência de diferentes graus e utilizando várias bases polinomiais. Para cada uma das formulações NNRPIM foi determinado o esquema de integração mais eficiente. Também os parâmetros de forma c e p da função de base radial multiquadrática foram analisados e otimizados.
As conclusões mais relevantes deste trabalho podem ser resumidas como:
- 1.
A formulação que apresenta resultados mais estáveis e precisos é a formulação NNRPIM v2p0.
- 2.
Para a formulação NNRPIM v2p0, os parâmetros de forma da função de base radial multiquadrática óptimos são c=0,0001 e p=0,9999.
- 3.
Determinou‐se que a eficiência da análise é maximizada se forem utilizados 2 pontos de integração em cada partição celular.
- 4.
Independentemente do tipo de problema em estudo, o NNRPIM apresenta consistentemente uma elevada taxa de convergência e uma elevada precisão. A solução final convergida é sempre muito próxima da solução analítica do problema considerado.
- 5.
O arranjo nodal não perturba significativamente a precisão ou a convergência do NNRPIM.
- 6.
A distribuição do campo de momentos e esforços transversos é suave ao longo de toda a viga. As perturbações nas fronteiras essenciais da viga devem‐se à menor conectividade nodal dos nós periféricos.
O NNRPIM provou ser um método alternativo válido para a análise de vigas laminadas.
Os autores agradecem o financiamento disponibilizado pelo Ministério da Ciência, Tecnologia e Ensino Superior – Fundação para a Ciência e Tecnologia (Portugal), através da bolsa SFRH/BPD/75072/2010 e pelo FEDER/FSE, através da bolsa PTDC/EME‐TME/098050/2008.