IEEJ TRANSACTIONS ON ELECTRICAL AND ELECTRONIC ENGINEERING IEEJ Trans 2020; 15: 157–158 Published online in Wiley Online Library (wileyonlinelibrary.com). DOI:10.1002/tee.23038 Letter Electromagnetic and Thermal Analysis of a Multilayer CFRP Panel Struck by Lightning with the FDTD Method Koki Ueno * , Student Member Nozomu Miki * , Student Member Yoshihiro Baba * a , Member Naoto Nagaoka * , Member Hiroyuki Tsubata ** , Non-member Takayuki Nishi ** , Non-member The transient distribution of heat generated in a 16-layer carbon fiber-reinforced plastic (CFRP) panel, in which a lightning current with a magnitude of 3 kA, a rise time of 6.4 μs, and a time to half-peak value of 69 μs is injected, has been computed. Each layer of the CFRP panel has a fiber axis direction of 45, 90, −45, or 0 ◦ and, therefore, has anisotropic electrical and thermal conductivities. This multilayer CFRP panel has been represented with a conductivity matrix consisting of off-diagonal elements in the finite-difference time domain (FDTD) simulation. The heat is computed, on the basis of FDTD-computed conduction current densities, using a discretized heat equation.  2019 Institute of Electrical Engineers of Japan. Published by John Wiley & Sons, Inc. Keywords: FDTD method; carbon fiber-reinforced plastics; Maxwell’s equations; heat equation Received 9 May 2019; Revised 1 July 2019 1. Introduction Carbon fiber-reinforced plastics (CFRP) have been recently used as the main material of aircraft bodies and wings. A CFRP panel is usually composed of many thin layers, each layer of which has a fiber axis direction of 45, 90, −45 or 0 ◦ . The electrical and thermal conductivities are high in the direction of fiber axis, and it is low in other directions. To the analysis of the distribution of lightning current in a CFRP panel, the finite-difference time domain (FDTD) method [1] has been applied [2,3]. In [3], a conductivity matrix with off-diagonal elements is used to represent a multilayer CFRP panel, while in [2], triangular-prism cells are used. The temperature of a CFRP-made aircraft might significantly increase locally when it is struck by lightning, which might cause partial deterioration and/or damage of the CFRP body and/or wing. In this letter, a procedure for computing heat [4] is incorporated in the FDTD simulation in the three-dimensional (3-D) Cartesian coordinate system, which uses a conductivity matrix with off- diagonal elements [3]. 2. Heat Computation Procedure The temperature T of a lossy material is obtained by solving the following heat equation: ∂ T ∂ t = α∇ 2 T + P d ρ m C m (1) where t is the time, α is the thermal diffusivity that is given by α = κ /(ρ m C m ), ρ m is the mass density, C m is the specific heat, κ is the thermal conductivity, and P d is the absorbed power per unit volume. In this letter, ρ m is set to 1.56 × 10 6 gm −3 , a Correspondence to: Yoshihiro Baba. E-mail: ybaba@mail.doshisha.ac.jp *Department of Electrical Engineering, Doshisha University, Kyotanabe, Kyoto 610-0394, Japan **SUBARU Corporation, Utsunomiya, Tochigi 320-8564, Japan C m is set to 0.944 J g −1 K −1 , the thermal conductivity κ in the fiber axis direction is set to 18 W m −1 K −1 , κ in the horizontally perpendicular direction is set to 1.9 W m −1 K −1 , and κ in the vertical direction is set to 0.95 W m −1 K −1 [5]. The electrical conductivities σ in the three above-mentioned directions are set to 33 000, 200, and 0.2 S m −1 , respectively, on the basis of the measurement. P d in a lossy material is expressed as follows: P d = J · E = (σ E ) · E (2) where J is the conduction current density vector, and E is the electric field vector. The absorbed power P d at the location (i , j , k ) is calculated by the following expression: P n+1 d (i , j , k ) = 1 t ′ ns +nt  m= ns [P m x (i , j , k ) + P m y (i , j , k ) + P m z (i , j , k )]t (3) where t is the field computation time increment; P x m (i ,j ,k ), P y m (i ,j ,k ), and P z m (i ,j ,k ) are instantaneous values at the time step m of the absorbed power in x , y , and z directions, respectively; n s is the first field computation time step number in a (long) thermal computation time increment t ′ ; and n t is the last field computation time step number in t ′ . The discretized equation of (1) is given below: T n+1 (i , j , k ) = αt ′ x 2 [T n (i − 1, j , k ) + T n (i + 1, j , k )] + αt ′ y 2 [T n (i , j − 1, k ) + T n (i , j + 1, k )] + αt ′ z 2 [T n (i , j , k − 1) + T n (i , j , k + 1)] +  1 − 2 αt ′ x 2 − 2 αt ′ y 2 − 2 αt ′ z 2  T n (i , j , k ) + P d n+1 (i , j , k ) · t ′ ρ m C m (4)  2019 Institute of Electrical Engineers of Japan. Published by John Wiley & Sons, Inc.