This is fleshed out from comments of Sam Lichtenstein; Some (or all) of what I have written is surely not the most elegant programming (at best), feel free to improve. One can do the following with MACAULAY2 (if you copy and paste this into the MACAULAY2 prompt it should work:)

R = GF(2)[x,y]/(x^3,y^3); m = ideal(x,y); ModSquare = ideal -> length(ideal/(ideal*ideal)); Generators = ideal -> length(ideal/(ideal*m)); I := ideal(x^2 + y^2); J := ann(I); [Generators(I),Generators(J),ModSquare(I),ModSquare(J)]

The function ModSquare applied to an ideal $I$ computes the length of $I/I^2$, and the function Generators computes the minimal number of generators of $I$. Similarly, ann computes the annihilator of $I$. From this we may compute that $I$ is principal, $J$ has two generators, and both $I/I^2$ and $J/J^2$ have length $4$. These functions only work for homogenous ideals. Following James Parson's suggestion, I also considered MAGMA (again with the restriction to affine algebras), and the following works for arbitrary ideals:

A:=AffineAlgebra<GF(2),x,y|x^3,y^3>; x:=A.1; y:=A.2; AssignNames(~A,["x","y"]); m:=ideal<A|x,y>; Minus := func<ideal | Dimension(quo<A|ideal>)>; Generators := func<ideal | Minus(ideal*m) - Minus(ideal)>; ModSquare := func<ideal | Minus(ideal*ideal) - Minus(ideal)>; ann := func<ideal | Annihilator(ideal)>; I:=ideal<A|x^2+y^2>; J:=ann(I); [Generators(I),Generators(J),ModSquare(I),ModSquare(J)];