How can I get the initial velocities and positions of the planets

NASA has an excellent site at http://ssd.jpl.nasa.gov/horizons.cgi#top. From there you can extract initial conditions in order to start your differential equation solver. At the above website you need to change from OBSERVER to VECTOR and then write in the planet you are interested in. The generated data contain the \( x \), \( y \) and \( z \) values as well as their corresponding velocities. The velocities are in units of AU per day. Alternatively they can be obtained in terms of km and km/s.

For the first simple system involving the Earth and the Sun, you could just initialize the position with say \( x=1 \) AU and \( y=0 \) AU.