Buried out some old pascal source :D About 14 years ago, I used it for displaying very simple 3d objects. xrot, yrot, zrot are to rotate points ([x,y,z] multiplied by rotation matrix). And I used a very simple 3d-to-2d-transformation, based on vanishing-point projection with the vanishing-point at the middle of the screen. As an example, there is a vertex array defined. You also have to add a trigons array.
const depth = 1500;
      deg = pi / 180;
      { some vertices for a dice :) }
      vertices:array[0..23] of real= (50, 50, 50,       { 0}
                                 -50, 50, 50,       { 1}
                                  50,-50, 50,       { 2}
                                 -50,-50, 50,       { 3}
                                  50, 50,-50,       { 4}
                                 -50, 50,-50,       { 5}
                                  50,-50,-50,       { 6}
                                 -50,-50,-50,       { 7}
                                  );
{ transform 3d coordinates to pixel coordinates }
procedure 3d_to_2d(x, y, z : real; var px, py : longint);
 var k:real;
begin
 k:=((depth shr 1)+z)/depth;
 px:=(getmaxx shr 1)+trunc(x*k);      { getmaxx is the width of the screen }
 py:=(getmaxy shr 1)+trunc(y*k);      { getmaxy is the height of the screen }
end;
{ rotate around the x axis by rx degrees }
procedure xrot(var x,y,z:real;rx:integer);
 var x1,y1,z1:real;
begin
 y1:=(y * cos(rx * deg))+(z* (sin(rx * deg)));
 z1:=(-y* sin(rx * deg))+(z* (cos(rx * deg)));
 y:=y1; z:=z1;
end;
{ rotate around the y axis by ry degrees }
procedure yrot(var x,y,z:real;ry:integer);
 var x1,y1,z1:real;
begin
 x1:=(x * cos(ry * deg))+(z*(sin(ry * deg)));
 z1:=(-x * sin(ry * deg))+(z*(cos(ry * deg)));
 x:=x1; z:=z1;
end;
{ rotate around the z axis by rz degrees }
procedure zrot(var x,y,z:real; rz:integer);
 var x1,y1,z1:real;
begin
 x1:=(x* cos(rz * deg))+(y*(sin(rz * deg)));
 y1:=(-x* sin(rz * deg))+(y*(cos(rz * deg)));
 x:=x1; y:=y1;
end;
For filled trigons, I used a friend's function, which draws the shape using horizontal lines (Hline(x,y,width, color)):
    TYPE pt=RECORD x,y:LongInt;END;  
PROCEDURE Tri(P:ARRAY OF pt;co:BYTE);
VAR q,w:INTEGER;
    S:pt;
    f12,f13,f23:LongInt;
    s1,s2:LongInt;
BEGIN
  IF p[0].y>p[2].y THEN BEGIN s:=p[0];p[0]:=p[2];p[2]:=s;END; { sort the points }
  IF p[0].y>p[1].y THEN BEGIN s:=p[0];p[0]:=p[1];p[1]:=s;END;
  IF p[1].y>p[2].y THEN BEGIN s:=p[1];p[1]:=p[2];p[2]:=s;END;
  q:=(p[0].y-p[1].y); { y distance between point 0 and 1 }
  IF q<>0 THEN f12:=LongInt((p[0].x-p[1].x) shl 6) DIV q ELSE f12:=0;
  q:=(p[0].y-p[2].y);
  IF q<>0 THEN f13:=LongInt((p[0].x-p[2].x) shl 6) DIV q ELSE f13:=0;
  q:=(p[1].y-p[2].y);
  IF q<>0 THEN f23:=LongInt((p[1].x-p[2].x) shl 6) DIV q ELSE f23:=0;
  s1:=p[0].x shl 6;s2:=s1;
  FOR q:=p[0].y TO p[1].y DO
  BEGIN
    Hline(s1 shr 6,s2 shr 6,q,co);
    s1:=s1+f12;
    s2:=s2+f13;
  END;
  s1:=p[2].x shl 6;s2:=s1;
  FOR q:=p[2].y DOWNTO p[1].y DO
  BEGIN
    Hline(s1 shr 6,s2 shr 6,q,co);
    s1:=s1-f23;
    s2:=s2-f13;
  END;
END;