Motivation: High-throughput experimental techniques have been producing more and more protein-protein interaction (PPI) data. PPI network alignment greatly benefits the understanding of evolutionary relationship among species, helps identify conserved sub-networks and provides extra information for functional annotations. Although a few methods have been developed for multiple PPI network alignment, the alignment quality is still far away from perfect and thus, new network alignment methods are needed. Result: In this paper, we present a novel method, denoted as ConvexAlign, for joint alignment of multiple PPI networks by convex optimization of a scoring function composed of sequence similarity, topological score and interaction conservation score. In contrast to existing methods that generate multiple alignments in a greedy or progressive manner, our convex method optimizes alignments globally and enforces consistency among all pairwise alignments, resulting in much better alignment quality. Tested on both synthetic and real data, our experimental results show that ConvexAlign outperforms several popular methods in producing functionally coherent alignments. ConvexAlign even has a larger advantage over the others in aligning real PPI networks. ConvexAlign also finds a few conserved complexes among 5 species which cannot be detected by the other methods.