Python, 193 191 bytes

F=lambda A,B,C:[p for i in[0,1,2]for j in range(i+1,4)if(z:=A[i]*B[j]-A[j]*B[i])if(p:=((C[i]*B[j]-C[j]*B[i])/z,(C[j]*A[i]-C[i]*A[j])/z))if min(A[k]*p[0]+B[k]*p[1]<C[k]+1e-9for k in range(4))]

Attempt This Online!

The arguments A, B, and C are the lists of coefficients in x, y, and the right-hand side respectively, i.e. A = [a, d, g, j], B = [b, e, h, k], and C = [c, f, i, l].

This is simply the exhaustive method of computing every intersection point and verifying if they satisfy every inequality. We first iterate over every pair of inequalities with for i in[0,1,2]for j in range(i+1,4), and we assign the variable z:=A[i]*B[j]-A[j]*B[i] using if. Since we assume no lines are parallel, z is always non-zero, and with some algebra we find that the intersection point is given by the coordinates p:=((C[i]*B[j]-C[j]*B[i])/z,(C[j]*A[i]-C[i]*A[j])/z), which we assign using if once more.

The last statement if min(A[k]*p[0]+B[k]*p[1]<C[k]+1e-9for k in range(4)) verifies that p is in the polygon by checking that it satisfies every inequality, which is the case if the "minimum" of the conditions is True. The +1e-9 is there to compensate for floating-point artefacts, since we also check the lines on which the point p lies. Note: the lack of a space between 1e-9 and for yields a SyntaxWarning, but the code still runs. Adding the space back in removes the warning at the cost of 1 byte.

This is my first post on here, please let me know if I forgot anything! This was fun, I learned a lot trying to chip away as much as I could :]

Edit: removed unnecessary brackets in the min.

Scala 3, 567 350 bytes

- 217 bytes thanks to @squareroot12621

Golfed version. Attempt This Online!

(A,B,C)=>{val P=scala.collection.mutable.ListBuffer[(Double,Double)]();val n=A.length;for(i<-0until n;j<-i+1until n){val a=A(i);val b=B(i);val c=C(i);val x=A(j);val y=B(j);val z=C(j);val d=a*y-x*b;if(d!=0.0){val q=(c*y-z*b)/d;val r=(z*a-c*x)/d;val p=(q,r);var v=true;for(k<-0until n if v)if(A(k)*p._1+B(k)*p._2>C(k)+1e-9)v=false;if(v)P+=p}};P.toList}

Ungolfed version. Attempt This Online!

  def findIntersectionPointsImproved(A: Seq[Double], B: Seq[Double], C: Seq[Double]): List[(Double, Double)] = {
    val intersectionPoints = ListBuffer[(Double, Double)]()
    val numLines = A.length

    for (i <- 0 until numLines) {
      for (j <- i + 1 until numLines) {

        val a_i = A(i)
        val b_i = B(i)
        val c_i = C(i)

        val a_j = A(j)
        val b_j = B(j)
        val c_j = C(j)

        val determinant = a_i * b_j - a_j * b_i

        if (determinant != 0.0) { 
          val x = (c_i * b_j - c_j * b_i) / determinant
          val y = (c_j * a_i - c_i * a_j) / determinant
          val point = (x, y)

          val tolerance = 1e-9

          var satisfiesAllConstraints = true
          breakable { 
            for (k <- 0 until numLines) {
              val constraintValue = A(k) * point._1 + B(k) * point._2
              if (constraintValue > C(k) + tolerance) {
                satisfiesAllConstraints = false
                break() 
              }
            }
          }

          if (satisfiesAllConstraints) {
            intersectionPoints += point
          }
        }

      }
    }
    intersectionPoints.toList
  }

Charcoal, 58 bytes

ＥΦＥΣＥθＥ…θκＥ³⁻קι⊕ν§λ⊖νקλ⊕ν§ι⊖ν∕ι§ι²⬤θ‹Σ×λιＸχ±χ⭆¹…ι²

Attempt This Online! Link is to verbose version of code. Takes input with negated constants as compared to the example i.e. ax + by + c <= 0 etc. Explanation:

ΣＥθＥ…θκ

Looping over every pair of lines...

Ｅ³⁻קι⊕ν§λ⊖νקλ⊕ν§ι⊖ν

... calculate their cross product in 3D space, ...

Ｅ...∕ι§ι²

... normalise to 2D points in homogenous coordinates, ...

Φ...⬤θ‹Σ×λιＸχ±χ

... filter on those inside the polygon (with a small fudge factor due to floating-point inaccuracy), ...

Ｅ...⭆¹…ι²

... and pretty-print the matching points.

Maple, 111 bytes

proc(s)uses PolyhedralSets;evalf(map2(eval,[x,y],Relations~(Vertices(PolyhedralSet(convert(s,rational))))))end;

Input is a set of the inequalities. Output is a list of coordinates as [x,y].

The task is just finding a 2D convex polytope's vertices from its H-representation. These operations are available in Maple's PolyhedralSets package. The package works in exact arithmetic, so convert(..,rational) is the first step. PolyhedralSet makes an internal data structure, Vertices finds the vertices, and Relations extracts them from the internal data structure. At this point we have a list of equations, e.g., [[y = -2379/878, x = -11952/2195],...]. The rest converts to a list of coordinates, e.g., [[-5.445102506, -2.709567198],...]. If [[y = -2.709567198, x = -5.445102506],...] were an acceptable output, the code could be shortened to 95 bytes:

proc(s)uses PolyhedralSets;(evalf@Relations)~(Vertices(PolyhedralSet(convert(s,rational))))end;

Lines going through the same point return that point, e.g., {x<=0,y<=0,x<=y} returns [[0., 0.]].