A new flux expansion nodal method is developed to solve the neutron diffusion equations in hexagonal-z geometry.The intra-nodal flux distribution is expanded in a series of analytic basis functions for each group.In order to improve the nodal coupling relations,a new type of nodal boundary conditions is proposed,which simultaneously requires the continuity of both the zero-and first-order moments of partial current across the nodal surfaces.The response matrix technique is used,which gives a fast-running scheme for the iterative solution of the nodal diffusion equations.Based on the proposed model,the code FEMHEX has been developed and tested against 2-D and 3-D benchmark problems for the VVER-type reactors.The numerical results show that FEMHEX can predict accurately the multiplication factor and nodal powers.