The effect of hawkmoth-like flexibility on the aerodynamic hovering performance of wings at a Reynolds number of 400 has been assessed by conducting fluid structure interaction simulations incorporating a finite difference based immersed boundary method coupled with a finite-element based structure solver. The stiffness distribution of a hawkmoth forewing was mapped onto three wing shapes (r¯1 = 0.43, 0.53, and 0.63) defined by the radius of the first moment of wing area each with aspect ratios, AR = 1.5, 2.96, 4.5, and 6.0 using elliptic mesh generation, the Jacobi method for iterations, and the concept of the barycentric coordinate system. The results show that there is a dominant chordwise deformation at AR = 1.5, and the wings also deform in the spanwise direction and their tips deviate from the horizontal stroke plane as AR increases. At AR = 1.5, 2.96, and 4.5, flexibility increases the mean lift (up to 39%, 18%, and 17.6%, respectively) for all wing shapes. At AR = 6.0, the r1¯ = 0.53 and 0.63 flexible wings give lesser lift than the rigid equivalents because of negative lift or small positive lift during the early stroke as the vortical structures remain on the bottom surface. This is attributed to the rapid pitch-down rotation, lesser stroke angular velocity than the rigid wing, and upward motion of the wingtip, away from the horizontal stroke plane. From the design perspective, the anisotropic flexible wings (except r1¯ = 0.53 and 0.63 with AR = 6.0) can be used in micro aerial vehicles for high lift requirements, such as for a high payload. Results here show that in nature, the hawkmoth wings with r1¯ and AR of 0.43-0.44 and 2.73-2.92, respectively, appear to have a combination of the shape, AR, and flexibility that optimizes power economy.