finite element – Implementing my own 2D FEM solver in C++

I’m trying to code a naive 2D FEM solver for my own understanding as a project in C++. Initially, I’m planning to start with heat equation with the plan of scaling it to other PDEs in the future. For meshes, I plan to use Gmsh and parse the mesh information.

It would be helpful, if someone can point out any resources for good code design specifically for fem solvers, or any documentation/tutorial that I can follow sequentially.
Since it’s a naive implementation, I’m interested in knowing about do’s and don’ts and not to overdo and get caught with complexities.